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SECTION  I 


INTRODUCTION 

Air  Force  sensor  and  weapon  requirements  for  Advanced  space-based  detec¬ 
tion,  surveillance,  and  defense  will  necessitate  the  construction  of  large 
space  structures  (LSS).  These  structures  will  be  very  large  and  lightweight, 
and  hence  extraordinarily  flexible.  They  will  most  probably  be  assembled  in 
space.  Their  generic  characteristics  would  include  the  interconnection  of 
several  rigid  bodies  with  flexible  elements,  with  the  rigid  bodies  also  sup¬ 
porting  flexible  appendages. 

Air  Force  mission  requirements  for  these  structures  will  impose  severe 
pointing  and  tracking  specifications.  Such  stringent  specifications  will 
necessitate  the  design  of  a  superior  control  system  that  would  involve  hun¬ 
dreds  of  sensors  and  actuators.  The  space  structure  dynamics  are  infinite 
dimensional,  and  design  models  may  involve  thousands  of  state  variables. 

Thus,  an  immense  amount  of  real-time  digital  computation  will  have  to  be 
carried  out  by  the  control  system.  This  massive  computation  cannot  be 
accomplished  by  a  centralized  control  computer.  Therefore  a  decentralized 
control  system  employing  distributed  computation  is  essential. 

The  need  for  decentralized  control  can  also  be  understood  by  the  fact 
that  these  very  large  structures  would  have  to  be  assembled  in  a  piecemeal 
manner.  A  partially  assembled  space  structure  will  need  some  real-time 
control  so  that  it  maintains  its  station,  performs  partial  functions,  and  does 


not  bend  during  further  assembly.  In  addition,  certain  LSS  (e.g.,  phased 
array  radars)  may  be  capable  of  functioning  with  only  a  partial  assembly.  The 
partial  assembly  will  require  a  complete  real-time  control  system  to  perform 
its  functions.  It  would  be  undesirable  to  completely  reprogram  in  space  a 
large  control  computer  every  time  an  addition  is  made  to  the  structure.  A 
more  sound  engineering  approach  is  to  augment  a  decentralized  control  system 
as  additional  elements  of  the  large  space  structure  (LSS)  are  brought 
together.  It  is  obvious  that  decentralized  and  hierarchical  control  archi¬ 
tectures  are  essential  if  these  envisioned  large  space  structures  are  to 
accomplish  their  mission. 

Another  factor  that-  significantly  influences  the  problem  of  control  sys¬ 
tem  design  for  large  space  structures  is  the  inherent  uncertainty  in  modeling 
the  dynamic  behavior  of  flexible  structures  in  space.  The  combination  of 
large  uncertainty,  large  dimensions,  and  severe  performance  specifications 
creates  a  difficult  control  system  design  problem  even  were  one  to  use  a 
centralized  control  architecture.  The  need  for  decentralized  and  hierarchical 
architectures  serves  to  transform  the  design  problem  into  one  which  has  re¬ 
ceived  minimal  attention  and  for  which  few  results  or  design  and  analysis 
procedures  have  been  obtained. 

The  primary  objective  of  this  project  has  been  to  initiate  research  to 
develop  a  sound  engineering  approach  to  the  design  of  decentralized  robust 
control  systems  for  large  space  structures.  The  engineering  procedure  to  be 
developed  must,  at  a  minimum,  address  the  following  issues: 

(a)  The  necessary  architecture  of  the  control  system  must  be  de¬ 
termined  from  the  system  structure  of  the  large  space  structure, 
the  performance  requirements  of  the  LSS,  hardware  implementation 
considerations,  and  any  sequential  assembly  requirements  for 
the  LSS; 
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(b)  The  procedure  must  be  able  to  readily  handle  the  large 
numbers  of  input,  output,  and  design  model  state  variables 
that  will  be  present  in  a  realistic  design  problem; 

(c)  The  procedure  must  be  able  to  incorporate  the  performance 
specifications  for  the  LSS  in  the  design  process; 

(d)  The  impact  of  dynamic  modeling  errors  due  to  the 
uncertainties  of  modeling  an  inherently  infinite  dimensional 
system  must  be  incorporated; 

(e)  The  impact  of  dynamic  modeling  errors  due  to  the  deliberate 
model  simplifications  that  will  be  made  to  accommodate  (a) 
and  (b)  must  be  incorporated. 

As  discussed  in  the  preceding  paragraphs,  the  required  control  system 
architecture  will  likely  require  both  decentralization  and  hierarchical 
control  layers  to  satisfy  the  requirements  listed  in  (a)  and  the  computational 
and  analytical  issues  raised  in  (b).  The  issues  raised  in  (c)  and  (d),  when 
taken  independently,  are  standard  issues  which  must  be  addressed  by  control 
system  designs  for  any  system.  The  issue  of  deliberate  modeling  error  is  also 
of  concern  in  control  system  design  for  many  systems.  However,  in  the  context 
of  LSS,  these  issues  assume  special  significance  when  combined  with  the  tight 
performance  specifications,  large  system  models,  and  non-classical  control 
architecture. 

Taken  together,  the  issues  raised  by  (a)-(e)  present  a  formidable  control 
system  design  and  analysis  problem  which  has  received  very  little  attention. 
The  combination  of  robustness  and  performance  with  a  decentralized  architec¬ 
ture  has  received  partial  consideration  only  in  [l]-[5].  References  [ 1 ]— [ 2 ) 
were  actually  concerned  with  the  development  of  decentralized  synthesis  pro¬ 
cedures,  but  the  results  can  also  be  used  to  analyze  the  robustness  of  weakly 
coupled  systems  (as  noted  in  [4]).  Reference  [3]  developed  a  design  procedure 
that  achieved  arbitrary  performance  specifications  using  a  decentralized 
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architecture  for  systems  that  had  no  bandwidth  limitations.  The  first  simul¬ 
taneous  consideration  of  performance,  robustness,  and  decentralized  archi¬ 
tectures  occurred  in  [4]-[5].  This  work  developed  a  preliminary  approach  to 
both  qualitative  and  quantitative  analysis  of  the  feedback  properties  of 
decentralized  systems.  Finally,  the  combination  of  robustness  analysis  with  a 
control  architecture  induced  by  a  system  time  scale  separation  structure  has 
been  considered  in  [6 ]— [7 ] . 

This  study  of  decentralized  robustness  procedures  for  large  space  struc¬ 
tures  attempts  to  organize  the  issues  discussed  in  the  preceding  paragraphs, 
identify  the  useful  analytical  tools  that  are  available  [l]-[6],  and  extend 
these  tools  where  possible.  Clearly  this  is  an  ambitious  undertaking.  Our 
initial  approach  limits  the  scope  of  the  research  to  the  study  of  a  specific 
simplified  problem  that  retains  representative  features  of  the  generic  de¬ 
centralized  design  problem  for  large  space  structures  with  sequential 
assembly.  The  limited  scope  allows  us  to  address  the  objectives  of  the  study 
in  a  well  defined  context,  thereby  allowing  the  development  of  insights  and 
understanding  into  the  Interactions  of  the  numerous  system  structure  and  con¬ 
trol  architecture  issues. 

Hence,  generic  system  structure  assumptions  (such  as  those  used  in  [1]- 
[7])  can  be  evaluated  for  reasonableness  and  applicability,  and  the  irrelevant 
structures  can  be  quickly  eliminated.  The  usefulness  of  the  resulting  analy¬ 
sis  tools  can  be  determined  in  terms  of  their  ability  to  produce  accurate  pre¬ 
dictions  of  system  performance  and  design  specifications. 

The  disadvantage  of  using  a  specific  system  as  a  basis  for  such  a  study 
is  the  difficulty  of  quantitatively  asserting  that  the  conclusions  of  the 


study  can  be  applied  beyond  the  specific  system. 
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This  study  minimizes  this 


difficulty  by  retaining  as  many  generic  features  of  the  large  space  structure 
problems  as  possible.  Although  the  quantitative  results  may  be  difficult  to 
generalize  precisely,  the  qualitative  assessments  of  available  analysis  design 
tools,  and  the  insights  into  the  interactions  between  system  structure,  con¬ 
trol  problem  structure,  and  control  system  architecture  will  be  valid  for  many 
problems. 

The  specific  problem  that  is  considered  is  the  sequential  assembly 
problem  for  two  subassemblies.  It  is  assumed  that  each  subassembly  consists 
of  a  highly  oscillatory  system  that  requires  vibrational  control  to  damp  the 
oscillations.  The  individual  control  systems  roust  perform  their  functions  on 
both  the  decoupled  subassemblies  and  on  the  interconnected  assembled  system. 

A  second  level  controller  using  a  reduced  set  of  controls  and  measurements  can 
be  implemented  after  the  interconnection  is  effected  to  augment  the  decentral¬ 
ized  controllers. 

This  particular  problem  was  chosen  for  its  ability  to  present  the  issues 
(a)-(e)  in  the  context  of  a  simplified  problem.  The  control  system  archi¬ 
tecture  (decentralized  control  with  a  second  level  coordinating  layer)  can  be 
(and  is)  deduced  from  the  system  structure,  performance  requirements  and  se¬ 
quential  assembly  requirements.  The  use  of  decentralized  control  directly 
addresses  the  problems  associated  with  large  numbers  of  variables.  The  fact 
that  each  subassembly  is  representative  of  a  LSS  implies  that  the  LSS  control 
issues  must  also  be  addressed  and  that  modeling  errors  that  result  from 
approximating  an  infinite  dimensional  system  must  be  considered.  Finally, 
both  modeling  errors  due  to  low  order  models  of  the  infinite  dimensional 
system  and  the  decentralized  design  procedures  are  incorporated. 
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Probably  the  most  significant  contributions  of  this  research  are  the 
identification  of  issues  and  problems  associated  with  the  decentralized  con¬ 
trol  of  large  space  structure  with  sequential  assembly,  and  the  development  of 
a  preliminary  design  methodology  for  such  systems.  At  a  more  technical  level, 
two  new  multivariable  robustness  results  have  been  developed.  Finally,  we 
have  also  identified  a  generic  error  structure  for  the  sequential  assembly 
problem. 

The  system  that  forms  the  basis  for  this  study  is  presented  in  Section  2. 
Both  the  subassemblies  and  the  interconnected  system  are  described  by  partial 
differential  equation  models.  Decomposition  results  for  such  systems  (pre¬ 
sented  in  Appendices  A  and  B)  are  used  to  analyze  the  structure  of  the  inter¬ 
connected  system  and  to  deduce  necessary  features  of  the  control  system 
architecture.  Low  order  design,  analysis,  and  truth  models  are  then  derived 
from  the  infinite  dimensional  models  using  low  order  finite  element 
approximations. 

Section  3  conducts  an  error  analysis  of  the  system  with  the  goal  of 
specifying  design  constraints  for  the  decentralized  control  systems.  Three 
generic  error  sources  are  identified,  and  quantitative  analyses  are  performed. 
The  results  of  these  analyses  are  used  to  specify  the  subassembly  control  sys¬ 
tems.  The  control  system  designs  are  presented  in  Section  4  and  the  validity 
of  the  analyses  of  Section  3  is  verified.  Section  5  summarizes  the  contri¬ 
butions  and  conclusions  obtained  from  this  study. 


SECTION  2 


SEQUENTIAL  ASSEMBLY  MODEL 

2.1  INFINITE  DIMENSIONAL  MODELS 

The  purpose  of  this  section  is  to  present  a  model  that  will  be  used  to 
study  the  applicability  of  the  various  decentralized  robustness  analysis  and 
design  tools  [l]-[7]  in  the  context  of  the  sequential  assembly  of  large  space 
structures.  The  model  has  been  developed  to  be  as  simple  as  possible  while 
retaining  the  basic  features  of  the  sequential  assembly  problem. 

A  number  qualitative  features  of  large  space  structures  and  the  sequen¬ 
tial  assembly  problem  distinguish  the  control  problems  for  such  systems  from 
other  control  problems.  Large  space  systems  are  inherently  infinite  dimen¬ 
sional  with  closely  spaced,  lightly  damped  modes.  This  feature  implies  that 
any  design  based  on  a  finite  dimensional  design  model  must  be  able  to  cope 
with  a  (possibly  large)  number  of  lightly  damped  unmodeled  modes  in  the  roll¬ 
off  region  of  the  design.  The  controls  and  measurements  for  the  system  are 
typically  "point"  controls  and  measurements,  with  significantly  fewer  controls 
than  states.  The  sequential  assembly  problem  imposes  a  decentralized  control 
architecture  on  the  overall  assembly.  Each  of  the  decentralized  controllers 
must  function  with  its  separate  subsystem  and  must  be  designed  to  maintain 
stability  when  connected  to  the  overall  assembly.  Finally,  a  second  level 
control  architecture  can  be  used  to  minimize  the  performance  degradation 
caused  by  the  interconnection  of  the  subsystem  assemblies. 
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The  model  that  will  be  used  to  represent  these  features  is  based  on  the 
physical  system  depicted  in  Fig.  2-1.  Each  subsystem  consists  of  three  unit 
masses  connected  by  two  uniform  beams.  The  two  subsystem  controls  are  the  , 
forces  on  the  middle  and  one  end  mass,  while  the  two  subsystem  measurements 
are  the  separations  of  the  masses.  The  generic  goal  of  the  subsystem  control¬ 
lers  is  to  maintain  the  nominal  separation  of  the  masses  as  nearly  as  pos¬ 
sible.  The  primary  contribution  to  errors  in  these  separations  are  longi¬ 
tudinal  vibrations  of  the  systems.  Assuming  that  the  masses  can  be  modeled  as 
point  masses  on  a  uniform  beam  of  length  L®,  the  longitudinal  position  of  an 
incremented  element  of  a  subsystem  can  be  described  by  the  hyperbolic  partial 
differential  equation 

32  3  32 

— “  z(x.t)  +  b  —  z(x,t)  -  c2  — -  z(x,t)  -  0 

3t2  3t  1  3x2  (2-1) 

where 


z(x,t)  *  longitudinal  position  of  the  point  x  of  the 
beam  at  time  t 

b  *  coefficient  of  restoring  force 
Cj  ■  velocity  of  wave  propagation 

The  boundary  condition  for  the  first  subsystems  are 


—  z(*,t) 
3x 


x-0 


3 

3x 


7“  z(*,t)  |l® 


32 

-  m  —  z(0,t)  +  un 
3t2 


t)  +  ui2 


(2-2) 


(2-3) 
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Figure  2-1.  Physical  Basis  for  Sequential  Assembly  Model 


32 

-  m  -  z(Lm,t) 

3t2 


while  the  boundary  conditions  for  the  second  subsystem  ares 

a  I  32 


—  z( *,t)  »  -m — -  z(0,t) 

3x  x-0  3t2 


—  z( •  ,t )  Lo  32  Lm 

3x  x* — '  -  -  m  -  z( — , t)  +  u2l 

2  3t2  2 


—  z(«,t) 
3x 


■  -  m -  z(L  ,t)  +  u22 


3x  lx-Lm  3t2 

where  m  is  the  mass  of  the  physical  masses  M1-M6. 


(2-4) 


(2-5) 


(2-6) 


(2-7) 


The  measurements  of  the  subsystems  are 


yn(t)  =  zi(— ,t) 

-  zi(0,t) 

(2-8) 

yi2(t)  -  zi(Lm,t) 

Lm 

-  zi(-.t) 

(2-9) 

Lm 

Y2l(t)  =  Z2(~»t) 

-  Z2(0,t) 

(2-10) 

y22(t)  -  z2(Lm,t) 

Lm 

-  *2Cpt) 

(2-11) 

For  each  subsystem,  the  parameters  b,  c,  Lm  and  m  were  chosen  as: 

b  =  .06 

(2-12) 

m  =  1 

(2-13) 

Lm  =  2 

(2-14) 

ci  =  1 

(2-15) 

The  subsystems  described  in  the  preceding  paragraph  will  define  the  sub¬ 
assemblies  for  the  sequential  assembly  problem.  The  connection  of  the  sub- 
assemblies  (with  length  L)  will  be  made  using  a  beam  of  length  Lm/2  whose 
longitudinal  wave  propagation  velocity  is  c2«  The  second  level  input  that 
will  be  used  to  maintain  the  separation  between  the  subassemblies  is  the 
differential  force  applied  between  the  end  masses  of  the  two  subassemblies, 
while  the  second  level  measurement  is  the  separation  of  those  masses.  The 
overall  system  is  described  by  the  non-uniform  hyperbolic  partial  differential 
equation 


32  3  o  32 

—  z(x,t)  +  b  —  z(x,t)  -  cz(x)  —  z(x,t)  =  0 

3t2  3t  3x2 


(2-16) 
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where 


3  ® 

—  z(0,t)  ““m  —  z(0,t)  +  ui  i 
3x  3t2 


3  L 

32  L 

—  z(-,t)  = 
3x  2 

”  n:  -  •  z(-,t)  +  ui? 

3t2  2 

(2-18) 

3 

32 

—  z(L,t)  = 
3x 

-  m  —  z(L,t)  + 

3t2 

(2-19) 

3  3 

32  3 

-  z^L.t)  - 

“  ®  — -  z(-L,t)  -  U3 

3t2  2 

(2-20) 

3 

32 

—  z(2L,t)  = 
3x 

~  m  —  z(2L,t)  +  u2i 

3t2 

(2-21) 

3  5 

32  5 

-  2(-L,t)  - 

“  m  *—  z(-L,t)  +  U22 
dt2  2 

(2-22) 

L 

yn(t) 

“  z(-,t)  -  z(0,t) 

(2-23) 

L 

yi2(t) 

*  z(L,t)  -  z(-,t) 

2 

(2-24) 

3 

y2i(t) 

-  z(2L,t)  -  z(-L,t) 

2 

(2-25) 

5 

y22<t) 

■  z(-h,t)  -  z(2L, t) 

(2-26) 

(  C1 

3  5 

0  <  x  <  L;  -L  <  x  <  -L 

2  2 

c(x)  - 

I 

1 

c2 

i 

3 

L  <  x  <  -  L 

2 

(2-27) 
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The  design  and  truth  models  will  be  derived  from  the  infinite  dimensional 
models  (2-1),  (2-16)  using  low  order  finite  element  methods.  The  use  of  low 
order  finite  element  models  allows  the  retention  of  most  of  the  features  of 
the  sequential  assembly  problem  and  large  space  structures.  However,  models 
derived  from  low  order  finite  element  approximations  are  typically  very  in¬ 
accurate  and  vary  strongly  with  the  order  of  the  approximation.  To  minimize 
the  effects  of  this  unrealistic  model  variation,  the  length  parameter  L  of  the 
assembled  subsystem  (2-16)  was  adjusted  to  yield  the  same  fundamental  mode  as 
the  separated  subsystem  (2-1)  when  C2  is  taken  to  be  zero.  Thus  the  subsystem 
models  generated  by  (2-16)  when  the  coupling  beam  is  not  present  will  be 
the  same  as  the  subsystem  models  generated  by  (2-1),  and  any  errors  between 
design  and  truth  models  can  be  attributed  to  physical  sources.  The  values  of 
L  that  are  used  will  depend  on  the  order  of  the  finite  element  approximation. 

2.2  STRUCTURAL  ANALYSIS  OF  THE  ASSEMBLED  SYSTEM 

The  longitudinal  wave  velocity  C2  of  the  coupling  beam  will  determine  the 
system  structural  properties  of  the  assembled  system  and,  consequently,  will 
affect  the  decentralized  and  hierarchical  structures  that  can  be  utilized. 
Three  cases  for  the  value  of  C2  will  be  qualitatively  examined  in  this  sub¬ 
section  (see  Fig.  2-2): 


CASE 

A: 

ci 

»  C2 

(2-28) 

CASE 

B: 

ci 

«  C2 

(2-29) 

CASE 

C: 

ci 

■  C2 

(2-30) 

In  both  cases  A  and  B,  the  asymptotic  eigenanalysis  of  Appendices  A  and  B 
and  the  analyses  of  [8]-[9]  apply.  The  analysis  for  case  A  leads  to  signifi¬ 
cant  conclusions  regarding  the  control  architectures  that  can  be  used  on  the 
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system.  Cases  B  and  C  are  more  difficult  to  address  precisely,  but  generic 
conclusions  can  be  drawn  from  similar  system  structures. 

Appendices  A  and  B  address  the  problem  of  spectral  analysis  and  solution 
approximation  for  a  class  of  stiff  systems.  This  class  is  defined  precisely 
in  the  appendices,  and  includes  hyperbolic  systems  of  the  form  of  (2-16)  with 
either  (2-28)-(2-29).  The  basic  conclusion  of  these  appendices  (and  the  re¬ 
lated  works  [9]— [ 10] )  applied  to  (2-16)  is  that,  for  the  purposes  of  analysis 
and  control,  (2-16)  can  be  decomposed  into  two  independent  analysis  or  control 
problems  (termed  the  stiff  and  normal  problems).  One  problem  is  based  on 
Infinite  dimensional  models  of  the  two  subassemblies,  while  the  other  problem 
approximates  (2-16)  by  a  single  beam  connecting  two  equivalent  masses  (i.e., 
the  subassemblies  are  approximated  by  equivalent  point  masses).  The  parti¬ 
cular  case,  A  or  B,  determines  the  forms  of  the  two  problems  and  hence  the 
control  architectures  that  will  be  effective. 

In  case  A,  (see  Fig.  2-2a),  the  interconnection  between  the  subassemblies 
is  much  less  stiff  than  the  subassemblies  themselves.  As  noted  in  the  pre¬ 
ceding  paragraph,  the  overall  problem  decomposes  into  two  Infinite  dimensional 
problems.  The  stiff  problem  approximates  (2-16)  as  two  independent  systems, 
each  of  which  is  modeled  by  (2—1)— (2— 11) .  Thus  any  controller  design  based  on 
(2—1)— (2— 11)  will  face  only  a  minimal  amount  of  error  when  applied  to  the 
coupled  system.  The  normal  problem  approximates  each  subassembly  by  an 
equivalent  point  mass  equal  to  the  total  subassembly  mass.  These  equivalent 
point  masses  are  connected  by  a  beam  whose  longitudinal  wave  velocity  is  C2« 
Thus,  this  structure  imposes  a  natural  two  level  decomposition  of  (2-16)  into 
two  decentralized  subassembly  controllers  and  one  coordinating  controller  that 
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regulates  the  subassembly  separation.  This  case  (2-28)  and  the  corresponding 
decentralized/hierarchical  architecture  will  be  investigated  in  detail  in  the 
remainder  of  this  report. 

Case  B  (sec  Fig.  2-2b)  also  decomposes  into  two  problems.  The  stiff 
problem  in  this  case  is  represented  by  two  equivalent  point  masses  connected 
by  a  beam  whose  longitudinal  wave  velocity  is  C2*  However,  the  equivalent 
masses  are  taken  as  the  end  masses  of  the  subassembly.  The  normal  problem  is 
based  on  the  subassemblies  but  does  not  decompose  as  in  Case  A.  Rather,  the 
model  assumes  the  form  of  (2-1)  with  5  point  masses  rather  than  3.  The  middle 
mass  represents  the  combined  mass  of  the  end  masses  of  the  subassembly.  Thus, 
the  normal  problem  is  an. almost  uniform  system,  with  the  only  non-uniformity 
being  the  doubled  center  mass. 

Consequently,  Cases  B  and  C  result  in  basically  the  same  control  prob¬ 
lems:  a  uniform  system  that  must  be  controlled  by  a  decentralized  control 
architecture.  The  connection  of  two  uniform  subassemblies  of  the  form  of  (2- 
1)  into  a  uniform  system  of  the  form  (2-16)  significantly  modifies  the  sub- 
assembly  transfer  function  (in  particular,  the  fundamental  modes  and  all 
harmonics  significantly  affected)  and  introduces  a  strong  coupling  between  the 
subsystems.  Any  significant  control  effects  on  the  subsystems  are  almost  sure 

to  push  one  or  more  of  the  lightly  damped  modes  into  the  right  half  plane. 

\ 

Hence  a  decentralized/hierarchical  structure  as  outlined  in  the  Introduction 
and  in  Subsection  2.1  is  not  feasible. 

Other  combinations  of  partially  decentralized  and  hierarchical  architec¬ 
tures  may  prove  useful  for  such  structures.  For  example,  the  use  of  over¬ 
lapping  local  controllers  has  proven  successful  in  controlling  lightly  damped 
systems  such  as  strings  of  moving  vehicles  [10]  and  freeway  systems  [11]. 


This  architecture  would  require  that  a  second  level  controller  which  overlaps 
the  subassembly  controllers  be  in  place  and  active  before  the  subassemblies 
are  connected.  Another  promising  approach  would  be  to  use  the  local  subas¬ 
sembly  controllers  to  induce  a  system  structure  that  is  amenable  to  decompo¬ 
sition  on  the  assembled  structure.  Finally,  additional  mechanisms  for  natural 
decompositions  of  the  assembled  system  exist  (e.g.,  large  end  masses  will  im¬ 
pose  a  decomposition  similar  to  Case  A).  These  mechanisms  and  their  effects 
can  be  identified  using  the  results  of  Appendices  A  and  B.  Each  of  these 
approaches  to  Cases  B  and  C  warrants  further  investigation. 

2.3  FINITE  ELEMENT  APPROXIMATIONS  FOR  DESIGN  AND  TRUTH  MODELS 

Low  order  finite  element  approximations  for  systems  (2-1)  and  (2-16)  were 
used  to  develop  design,  analysis,  and  truth  models  for  use  in  the  development 
and  illustration  of  the  decentralized/sequential  assembly  design  procedure. 

Due  to  the  strong  dependence  of  the  system  description  on  the  number  of  ele¬ 
ments  used  in  low  order  approximations,  some  of  the  parameters  of  the  design 
model  and  the  coupled  system  model  were  adjusted  to  reduce  the  difference 
between  the  fundamental  modes  of  the  models.  These  adjustments  affect  only 
the  errors  due  to  the  low  order  approximations,  and  not  the  qualitative 
properties  of  the  true  error  sources.  The  specific  adjustments  that  were  made 
are  indicated  in  the  description  of  the  individual  models. 

The  finite  element  method  that  was  used  employs  overlapping  triangular 
elements  (see  Fig.  2-3).  The  elements  on  each  end  of  the  beam  being  approxi¬ 
mated  are  right  triangles  whose  base  is  half  the  base  of  the  interior  ele¬ 
ments.  The  interior  elements  are  isosceles  triangles.  Each  interior  element 
overlaps  half  the  preceding  element  and  half  the  succeeding  element.  The 


Figure  2-3.  Example  of  5  Element  Approximation 


elements  are  normalized  with  respect  to  the  longitudinal  wave  velocity,  i.e., 

L 

/  c(x)  Vi  (x)2  dx  ■  1 
0 

Data  for  each  of  the  models  can  be  found  in  Appendix  C. 

The  truth  models  for  the  unassembled  systems  were  derived  from  (2-1)- 
(2-15)  using  a  5  element  approximation.  The  resulting  10th  order  system 
was  reduced  to  an  8th  order  system  by  eliminating  the  unobservable  center  of 
mass  (CM)  position  and  velocity  states.  The  resulting  truth  model  open  loop 
modes  were: 

-.03  ±  4 . 627 j 
-.03  ±  4.62>j 
-.03  ±  1 .485 j 

-.03  ±  .99j  (2-31) 


The  singular  values  of  the  truth  subsystem  model  transfer  function  are  shown 
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Figure  2-4.  Singular  Values  of  the  5  Element  Subsystem 
Truth  Model 

The  design  models  for  the  subassemblies  were  obtained  from  (2-1 )-( 2-11) 
using  a  3  element  approximation.  The  longitudinal  wave  velocity  was  adjusted 
to  align  the  fundamental  mode  of  the  design  model  with  the  truth  model.  The 
Oth  order  unobservable  modes  were  again  eliminated.  The  eigenvalues  of  the 
resulting  4th  order  design  systems  were: 

-.03  ±  1 . 4685 j 

-.03  ±  . 99 j  (2-32) 

The  singular  values  of  the  design  subsystem  model  transfer  function  are  shown 
in  Fig.  2-5. 

An  analysis  model  was  derived  from  the  coupled  system  (2— 12 )— (2— 27 )  using 
a  6  element  approximation  for  the  purpose  of  analyzing  the  perturbational  ef¬ 
fects  of  coupling  the  two  subassemblies.  The  subassembly  lengths  were 
adjusted  to  be 
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Figure  2-5.  Singular  Values  of  the  3  Element  Subsystem  Design 
Model 

tT 

S 

to  obtain  consistency  between  the  decoupled  truth  model  and  the  subassembly 
truth  model.  The  longitudinal  wave  velocity  of  the  coupling  beam  was  taken  to 
be 


c2  "  •! 


(2-34) 


The  two  eigenvalues  corresponding  to  the  Oth  mode  were  found  to  be  unobserv¬ 
able  and  were  removed.  The  resulting  10th  order  system  had  the  following 

% 

eigenvalues : 


-.03  ±  1.468 j 
-.03  ±  .9913 j 
-.03  ±  1.4669J 
-.03  ±  .9886j 
-.03  ±  .0542j 


subsystem  1  modes 


subsystem  2  modes 
coupling  mode 


(2-35) 
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Finally,  the  truth  uiodel  for  the  assembled  system  was  derived  from  (2- 
12)— (2— 27 )  and  (2-33)-(2-34)  using  an  11  element  approximation.  The  two 
eigenvalues  corresponding  to  the  Oth  order  mode  were  found  to  be  unobservable 
and  were  removed.  The  resulting  20th  order  system  had  the  following  eigen¬ 
values  : 


-.03 

+ 

4.2572j 

|  coupling  mode 

-.03 

± 

4.7  586 j  | 

1 

?  subsystem  1  unmodeled 

-.03 

4 

4. 5234 j 

j  dynamics 

-.03 

4 .7586 J  ( 

subsystem  2  unmodeled 

-.03 

4 

4.5234j 

dynamics 

-.03 

± 

1.3995J  j 
1 

1 

>  subsystem  1  modes 

-.03 

4 

.9397 j  ] 

1 

-.03 

± 

1.3995J  J 

1 

•  subsystem  2  modes 

-.03 

4 

•9397j  j 

1 

-.03 

± 

.0554 j  i 

coupling  mode 

(2-36) 
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SECTION  3 


ERROR/ STRUCTURAL  ANALYSIS  AND 
CONTROL  OBJECTIVES 


3.1  METHODOLOGY 

The  purpose  of  this  section  is  to  describe  a  methodology  for  the  design 

of  decentralized  control  systems  for  large  space  structures  that  are  to  be 

* 

assembled  sequentially,  and  to  develop  the  necessary  analytical  tools  to  im¬ 
plement  this  methodology.  The  framework  that  will  be  used  for  the  methodology 
is  common  to  all  control  system  design  problems  (see  Fig.  3-1).  The  first 
step  is  to  identify  all  performance  criteria  and  error  sources.  These  are 
then  mapped  (using  appropriate  analytical  tools)  into  constraints  on  the  con¬ 
trol  system  design.  The  control  system  must  be  designed  to  satisfy  the  re¬ 
sulting  constraints,  and  its  performance  is  verified  via  analysis  and  simula¬ 
tion.  Each  transition  between  steps  of  this  procedure  usually  requires 
iterations  to  guarantee  that  consistent  specifications  for  the  next  step  are 
achieved. 

Our  implementation  of  the  methodology  will  be  based  on  loop  shaping  syn¬ 
thesis  techniques  [12]— [13]  that  utilize  singular  value  analyses.  A  decen¬ 
tralized  design  will  be  developed  by  requiring  each  subsystem  to  use  a  local 
design  model.  The  subsystem  design  problems  will  then  be  determined  by  speci¬ 
fying  control  system  design  constraints  in  terms  of  bounds  on  the  singular 
values  of  the  loop  transfer  functions  of  the  closed  loop  subsystem  design 
models. 
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Figure  3-1.  Generic  Methodology  for  Control  System  Design 

Once  the  subsystem  design  problems  are  formulated  in  this  manner ,  the 
synthesis  of  a  subsystem  controller  that  satisfies  the  singular  value  con¬ 
straints  is  a  standard  problem  that  has  been  well  studied  (c.f.  [12]— [13]). 

The  unique  aspects  of  this  approach  for  decentralized  LSS  with  sequential 
assembly  are  the  type  of  error  sources  that  must  be  addressed,  and  the  need  to 
develop  singular  value  specifications  that  are  not  so  conservative  that  the 
resulting  synthesis  problem  is  unsolvable. 

The  main  emphasis  of  this  section  will  be  the  development  and  evaluation 
of  analysis  tools  that  map  the  unique  large  space  structure/ sequential  as¬ 
sembly  decentralized  control  problem  error  sources  into  specifications  on  the 
singular  values  of  the  subsystem  design  problems.  There  are  three  dominant 
error  sources  for  this  problem:  the  usual  system  modeling  errors,  a  modifica¬ 
tion  of  the  subsystem  transfer  function  model  due  to  the  neglected  coupling, 
and  interactions  between  the  subsystems  due  to  the  neglected  couplings.  An 
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appropriate  analysis  tool  is  developed  or  identified  for  each  error  source  to 
guarantee  that  the  resulting  control  system  design  problems  are  feasible. 

The  appropriateness  of  an  analysis  tool  for  a  particular  error  source 
depends  largely  on  the  presence  of  an  underlying  system  structure  for  the  LSS 
and  error  source.  The  presence  of  such  a  structure  presents  an  opportunity  to 
exploit  that  structure  to  reduce  the  conservativeness  of  the  corresponding 
design  specification.  In  fact,  large  error  sources  can  make  exploitation  of 
structure  essential.  One  of  the  significant  contributions  of  this  section  is 
to  identify  a  generic  structure  for  a  dominant  error  source  and  to  develop  an 
analysis  tool  that  exploits  this  structure. 

3.2  SIMPLIFIED  CONTROL  OBJECTIVES  AND  ERROR  SOURCES 

As  noted  in  the  preceding  subsection,  the  purpose  of  this  section  is  to 
develop  or  identify  analytical  tools  that  can  be  used  to  implement  the  decen¬ 
tralized  design  methodology  illustrated  in  Fig.  3-1.  To  concentrate  on  this 
objective,  the  methodology  will  be  applied  to  system  models  derived  from  the 
system  described  in  Section  2  with  simplified  performance  goals  and  error 
analyses.  The  simplifications  serve  to  concentrate  our  effort  on  the  method¬ 
ology  development  and  analysis  by  minimizing  the  conceptually  unimportant  de¬ 
tails  which  result  from  the  full  problem.  Meanwhile,  the  use  of  a  system  such 
as  that  of  Section  2  allows  us  to  retain  the  generic  structure  of  the  models 
and  error  sources. 

The  primary  objective  of  a  feedback  control  system  is  to  maximize  perfor¬ 
mance  in  an  appropriate  sense  within  the  constraints  Imposed  by  the  actuator 
and  sensor  limitations  and  plant  modeling  errors.  In  general,  the  performance 
goals  will  be  quantified  in  a  number  of  ways.  Examples  include  reduction  in 


the  system  response  due  to  classes  of  external  disturbances,  percent  accuracy 
in  following  classes  of  commands,  and  maintenance  of  such  goals  in  the 
presence  of  parameter  variations.  All  these  are  compatible  goals  and  serve 
primarily  to  quantify  the  system  properties  that  are  required.  Since  our  goal 
is  to  focus  on  a  systematic  procedure  for  handling  the  multiple  sources  of 
modeling  error  in  the  sequential  assembly  problem,  we  will  use  a  simplified 
performance  goal:  maximize  the  system  bandwidth  while  maintaining  stability 
and  zero  steady  state  errors  in  the  mass  separations.  We  will  view  this  goal 
in  a  loose  sense  in  that  we  will  also  attempt  to  evaluate  the  conservativeness 
of  the  specifications  that  result  from  the  error  analysis.  The  latter  will  be 
achieved  by  synthesizing  a  control  system  that  "nearly"  violates  the  design 
specifications. 

The  error  sources  that  will  be  considered  are  those  that  can  be  directly 
attributed  to  the  special  features  of  the  LSS/sequential  assembly  decentral¬ 
ized  control  problem.  These  error  sources  can  be  represented  by  differences 
between  the  design,  analysis  and  truth  models.  The  errors  can  be  associated 
with  three  sources,  each  of  which  must  be  addressed  differently  due  to  the 
structure  of  the  error  source. 

The  first  error  source  is  the  mismatch  between  the  local  design  model  for 
a  subassembly  and  the  actual  subassembly  dynamics.  For  purposes  of  simpli¬ 
city,  this  error  source  will  be  specified  as  the  error  between  the  subassembly 
design  and  truth  models  described  in  Section  2.  In  practice,  this  error 
source  would  have  to  be  quantified  to  account  for  all  possible  subsystem 
parameter  variations  and  for  all  the  high  frequency  unmodeled  dynamics.  How¬ 
ever,  the  general  character  of  this  source  can  be  represented  by  our  approach 


as  long  as  the  limitations  of  the  quantitative  error  analysis  based  on  this 
simplified  model  are  recognized. 

The  connection  of  the  subassemblies  leads  to  the  remaining  two  error 
sources.  The  interconnection  causes  the  subsystem  transfer  functions  (i.e., 
the  transfer  function  from  the  inputs  to  outputs  of  a  single  subsystem)  to 
change,  and  introduces  a  coupling  transfer  function  between  subsystems.  Al¬ 
though  the  control  system  design  is  based  on  the  subsystem  design  models,  it 
must  also  stabilize  the  interconnected  system.  This  means  that  the  nominal 
design  must  be  robust  to  the  cwo  error  sources  (subsystem  modification  and 
neglected  coupling  interactions). 

The  error  source  due  to  subsystem  modification  will  be  specified  as  the 
difference  between  the  subassembly  design  model  and  the  subassembly  transfer 
function  that  results  from  the  analysis  model.  The  general  approach  to 
characterizing  this  error  source  would  use  an  analysis  model  of  the  coupled 
system  that  is  at  least  as  detailed  as  the  design  models  for  the  subsystem. 

The  structure  and  magnitude  of  the  error  source  would  be  determined  from  the 
difference  between  the  analysis  subsystem  models  and  the  design  models.  The 
errors  would  then  be  characterized  in  terms  of  specifications  on  the  singular 
values  and  singular  subspaces  of  the  closed  loop  design.  We  will  follow  this 
approach  with  the  only  difference  being  the  level  of  detail  of  the  design  and 
analysis  models.  It  should  be  emphasized  that  the  analysis  model  is  only  used 
to  characterize  the  error  source  for  robustness  analysis,  and  impacts  the  con¬ 
trol  system  synthesis  only  through  this  characterization.  Hence  the  addition¬ 
al  complexity  of  this  model  does  not  automatically  translate  into  comparable 
controller  complexity,  and  the  decentralized  control  architecture  will  be 


retained. 


The  third  error  source,  the  neglected  subsystem  interactions  of  the 
interconnected  system,  will  also  be  specified  by  the  intersubsystem  coupling 
transfer  functions  of  the  analysis  model.  Again,  one  would  generally  use  a 
more  detailed  analysis  model  to  characterize  the  structure  and  magnitude  of 
this  error  source.  However,  the  level  of  detail  of  the  analysis  model  does 
not  conceptually  effect  this  approach.  As  with  the  subsystem  modification 
error,  this  error  characterization  will  be  used  to  determine  specifications 
for  the  subsystem  designs.  Since  the  designs  will  still  be  based  on  the  un¬ 
coupled  design  models,  the  resulting  control  system  will  retain  the  desired 
decentralized  architecture. 

3.3  ROBUSTNESS  ANALYSIS.  TOOLS 

3.3.1  General  Robustness  Theorems 

The  purpose  of  this  subsection  is  to  present  and,  when  necessary,  develop 
the  analysis  tools  that  are  required  to  implement  the  methodology  outlined  in 
the  preceding  subsections.  The  results  are  presented  in  a  uniform  framework 
based  on  a  restatement  of  a  singular  value  robustness  result  by  Lehtomaki  et. 
al.  [1A]  and  a  new  result,  Theorem  3.2,  that  allows  the  combination  of 
different  robustness  tests.  A  series  of  tests  which  fit  this  category  and  are 
applicable  to  the  problem  described  in  Section  2  are  then  presented  in  Sub¬ 
sections  3. 3. 2-3. 3. A.  , 

Let  Po(s)  be  the  nominal  model  of  the  plant  transfer  function,  and  let 
Pj-( 8)  be  the  true  plant  transfer  function.  Let  A(s)  denote  the  input 
multiplicative  error  betweea  Pq(s)  and  Pt(s): 


He  will  assume  that  Pt(s)  is  only  known  to  be  in  a  convex  set  P.  The  class 
of  allowable  perturbations  V  is  then  defined  as  the  convex  set: 

V  =  (A(s) :  A(s)  =  P^l(s)  [Pt(s)  -  P0(s)],  Pt(s)  e  P}  (3-2) 

[ 

Note  that  in  practice  V  may  be  either  specified  directly  or  derived  from  the 
specification  of  P. 

Let  K(s)  denote  the  feedback  compensation.  Then,  the  loop  transfer 
function  of  the  nominal  feedback  system  at  the  plant  input  will  be  defined 
by: 

L0(s)  =  K(s)  P0(s)  (3-3) 

and  the  loop  transfer  function  of  the  true  feedback  system  will  be  defined 
I  by: 

Lt(s)  =  K(s)  Pt(s)  (3-4) 

It  will  be  assumed  that  the  feedback  compensator  is  square  and  invertible. 

Hence  the  relationship  (3-2)  between  the  error  class  V  and  plant  class  P  re¬ 
mains  valid. 

\  We  will  assume  that  the  nominal  plant  model  Po(s)  and  each  plant  Pt(s)  e  P 

satisfy  the  following  assumptions:. 

(Al) 

(A2) 

i 
i 

(A3) 

(A4) 
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?0(s)  has  a  finite  dimensional  realization. 

Pt(s)  is  a  matrix  of  functions,  each  analytic  in  the  closed 
right  half  plane  except  for  a  finite  number  of  poles. 

Pt(s)  and  ?o(s)  have  the  same  number  of  closed  right  half 
plane  poles. 

[I  +  Po(s)]“l  is  analytic  and  bounded  in  the  closed  right  half 
plane  (i«e.,  stable). 


(A5)  If  Pt(s)  has  a  pole  at  s=jmo,  then  Po(s)  has  a  pole  at  s=juo. 

(A6)  P“l(s)  exist  except  at  a  finite  number  of  points  in  the 
0 

complex  plane. 

Finally,  let  Dr  denote  the  Nyquist  contour  which  closes  in  the  right  half 
plane  with  a  circle  of  radius  R  and  indents  into  the  left  half  plane  around 
each  jo>  axix  pole  of  Lq(s)  with  a  circle  of  radius  1/R. 

With  this  background  we  can  state  the  following  fundamental  robustness 
theorem. 

Theorem  3. 1 :  Assume  that  (A1)-(A6)  hold  and  that  there  exists  an  Rg  such 
that  for  all  R  >  Ro,  for  each  s  e  Dr  and  for  each  A(s)  e  V 

a  [I  +  L-l(s)  +  aA(s)]  >0  0  <  a  <  1  (3-5) 

_  0 

Then  the  true  closed  loop  system  sensitivity  matrix  [I  +  Lt(s)]-1  and  transfer 
function  matrix  [I  +  Lt(s)]"^  ^(s)  are  analytic  and  bounded  in  the  closed 
right  half  plane  (i.e.,  the  closed  loop  system  is  stable). 

Proof:  This  theorem  is  a  slightly  generalized  restatement  of  Theorem 
2.2  in  [14].  First,  ([14],  Theorem  2.2)  can  be  extended  to  distributed  sys¬ 
tems  by  invoking  the  generalized  Nyquist  criterion  [15]  in  the  proof  of 
([14],  Theorem  2.2).  Assumptions  (A1)-(A5)  imply  that  assumptions  a)-c)  of 
([14],  Theorem  2.2)  are  satisfied.'  Assumption  (A2)  implies  that  Lt(s)  is  a 
member  of  the  class  of  distributed  systems  for  which  the  generalized  Nyquist 
criterion  [15]  holds  (i.e.,  Lt(s)  e  B(o);  see  [15]  for  details).  Finally, 
(3-5)  is  simply  a  restatement  of  assumption  d)  of  ([14],  Theorem  2.2),  as  the 


following  argument  shows. 


The  condition 


det  {I  +  (I-a)  Lq(s)  +  a  Lt(s)]  -  0  Lt(s)  e  K P  (3-6) 

holds  for  some  value  s  e  Dr,  some  Lt(s)  e  P  and  some  a  e  [0,1]  if  and  only  if 

det  [I  +  Lo(s)-1  +  e  A(s)]  =  0  (3-7) 

for  the  same  s  e  Dr,  a  e  [0,1]  and  for  some  A(s)  e  V,  But  condition  (3-7) 
holds  if  and  only  if* 

jj  [I  +  Lq(s)“1  +  a  A(s)]  =  0  (3-8) 

Hence  (3-3)  holds  if  and  only  if  d)  of  ([14],  Theorem  2,2)  holds,  ^ 

This  Theorem  (and  its  parent  [14])  are  the  basis  for  most  of  the  guaran¬ 
teed  robustness  tests  to  data.  The  results  for  multiplicative  perturbations, 
[12],  [6]  additive  perturbations,  [12],  divisive  perturbations  [14],  phase  in¬ 
formation  [16],  and  weak  coupling  [l]-[2]  can  all  be  shown  to  guarantee  that 
(3-5)  (or  a  minor  variant)  holds  for  each  frequency  on  the  Nyquist  contour  Dr. 
Two  of  these  results,  multiplicative  perturbations  and  weak  coupling,  will  be 
used  in  subsequent  subsections,  and  a  third  result  for  perturbations  with  a 
specific  structure  will  be  derived  directly  from  Theorem  3.1. 

*Note  that  the  use  of  the  smallest  singular  value  In  (3-8)  can  be  replaced  by 
any  inverse  norm  function,  i.e.,  any  function  t  defined  as: 

t[A]  =  II  A"*  I”1 

where  II  •  II  is  any  induced  norm. 
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A  difficulty  that  arises  When  using  a  number  of  diverse  stability  tests 
that  are  satisfied  only  on  part  of  the  Nyquist  contour  is  that  Theorem  3.1  and 
the  results  [6],  [12],  [14],  [16]  do  not  apply  directly.  The  following 
Theorem  allows  the  integration  of  tests  to  guarantee  that  Theorem  3.1  holds 
and  that  the  resulting  true  closed  loop  system  will  be  stable. 

Theorem  3.2;  Let  denote  a  test  that  is  true  if  and  only  if 

condition  (3-5)  is  satisfied  for  each  s  e  %c  Dr,  and  each  A(s)  e  P.  Assume 

that  assumptions  (A1)-(A6)  hold,  that  {Ti(%)}N  are  a  set  of  tests  that  are 

N  1=1 

all  true,  and  that  U  ^i  “  Dr.  Then  the  closed  loop  system  is  stable. 

i=l 

Proof:  Given  any  s  e  Dr,  there  is  a  set  flj  such  that  s  e  ftj.  Then  the 
assumption  that  Tj(ftj)  is  true  implies  that  (3-5)  holds.  Hence  Theorem  3.1 
implies  the  closed  loop  system  is  stable. 

■ 

The  power  of  this  theorem  lies  in  its  ability  to  combine  different  types 
of  robustness  tests  over  different  frequency  ranges  into  a  comprehensive 
analysis  which  guarantees  stability  of  the  closed  loop  system.  For  example, 
the  general  multiplicative  modeling  error  test  [12]  may  be  a  tight  constraint 
at  high  frequencies,  but  conservative  at  low  frequencies.  If  a  structured 
analysis  test  can  be  developed  for  the  low  frequency  error,  then  It  can  be 
combined  with  the  high  frequency  multiplicative  error  analysis  to  yield  a  com¬ 
plete  robustness  result. 

\ 

For  the  system  described  in  Section  2,  three  error  sources  dominate.  As 
will  be  seen  in  subsection  3.4,  the  local  modeling  error  is  predominantly  a 
high  frequency  error  for  which  the  unstructured  modeling  error  analysis  is 
appropriate.  The  subsystem  model  error  that  results  from  connecting  the  sub- 
assemblies  dominates  at  low  frequencies  as  subsection  3.5  demonstrates.  This 
error  requires  a  structural  error  analysis  to  avoid  being  overly  conservative. 
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Finally,  the  neglected  subsystem  Interaction  that  results  from  connecting  the 
subassemblies  is  examined  in  subsection  3.6.  This  error  source  is  largest  in 
the  mid-frequency  range,  and  can  be  analyzed  using  a  weak  coupling  analysis 
analogous  to  that  developed  by  Bennett  and  Baras  [i]. 

Each  of  these  three  robustness  tests  that  will  be  used  subsequently  in 
this  report  will  be  presented  in  the  remainder  of  this  section.  In  each  case, 
we  will  use  the  identity 

T0(s)  -  II  +  Lo(s)-!]-1  (3-9) 

3.3.2  Unstructured  Multiplicative  Error 

The  general  multiplicative  error  robustness  test  has  been  developed  in 
[12]  and  has  been  used  in  a  variety  of  contexts.  The  unstructured  multiplica¬ 
tive  error  Am(s)  e  Vn  characterized  by  the  bound  jo>) : 

K  -  Um(s):  o[A(jw)]~l  >  (3-10) 

Stated  in  the  form  of  Theorems  3.1  and  3.2  the  unstructured  multiplicative 
error  robustness  test  for  a  frequency  range  Qm  is: 

Theorem  3.3:  (Unstructured  Multiplicative  Error  Robustness  Test)  If 

\ 

o  [T0(ju)]  <  fm(ja>)  V  a)  e  (3-11) 

then  condition  (3-5)  holds  for  all  u>  e 

Proof:  Condition  (3-5)  can  be  written  as 


1  [I  +  Lq( ]  -  a  a[An(ju))]  >  0 


a  e  [0,1] 


(3-12) 


Condition  (3-1 1)  is  obtained  by  coving  the  second  term  to  the  right  side  of 
the  inequality,  maximizing  over  a  (a-1),  inverting  both  sides  of  the  in¬ 
equality,  and  using  (3-9). 


3.3.3  Rank  One  Errors 

The  subsystem  model  error  that  results  from  connecting  the  subassemblies 
is  highly  strucured.  The  principal  structural  features  of  the  error  source 
A8(8)  are  that  As(s)  has  rank  one  and  that  As(s)  is  dominant  over  low  frequen¬ 
cies  (see  subsection  3.5  for  a  more  complete  discussion).  The  first  feature 
implies  that  only  one  of  the  singular  values  of  As(s)  is  non-zero.  Hence  we 
can  write  the  singular  value  of  decomposition  of  As(s)  as  (where  the  explicit 
dependence  on  s  has  been  suppressed  for  notational  clarity): 


D  l  VH 


(3-13) 


where 


0}  0 

J>  0_ 


U  -  [ui  U2  3 

V  -  (vi  V2I 


(3-14) 


Quantities  that  will  be  useful  in  the  statement  and  proof  of  Theorem  3.4  are 

\ 

the  angle  <j>  between  the  left  and  right  singular  subspaces  corresponding  to  01 
and  the  relative  complex  angle  0  between  the  left  and  right  singular  vectors 
(see  [17]— [18]  for  details).  These  angles  are  defined  uniquely  whenever  uj  is 
not  orthogonal  to  vj  by: 


viH  ui  *  cos  $  e3®  ()<$<-,  -ir  <0<ir 


(3-15) 


kV£->\‘ 


ir 

when  ui  is  orthogonal  to  Vj,  0  is  undefined  and  ♦  =  -  . 

The  second  structural  feature  of  As(s)  implies  that  .the  nominal  closed 
loop  function  can  be  approximated  by  the  identity  matrix  over  the  dominant 
freqency  range  fl8  of  the  structural  error 

a  (To(j«)  -  I]  <  e  <a  e  fts  (3-16) 

with  e  <  1.  Theorem  3.4  is  a  guarantee  that  (3-5)  holds  for  error  sources 
with  this  structure  and  a  given  value  of  e  in  (3-16). 

Theorem  3.4:  (Rank  One  Error  Robustness  Test).  Assume  that  the  nominal 
closed  loop  transfer  function  To(s)  satisfies  (3-16).  Define 

e 

r  *  —  +  sup  [<Ji  sin  $] 

1-e  ajsJlg 

A8eP8  (3-17) 

If  for  each  i8c  D 

|l  +  oi  co«  ♦  e J ® |  >  r  V  a)  e  ft8  (3-18) 

then  condition  (3-5)  holds  for  each  w  e  fl8. 

Proof:  Let  w  e  fl8  and  Ds  e  V8.  Since  To(jm)  satisfies  (3-14),  To”^(ju>) 
satisfies 

% 

a  [T“l( jui)  -  I]  <  —  (3-19) 

“  0  1-e 

Condition  (3-5)  can  be  rewritten  as: 

o  U'r1(ju>)  -  I]  +  I  +  aAs( jo))}  >  0  Voe[0,l]  (3-20) 
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Inequality  (3*20)  will  be  satlsifed  if 


£  [I  +  a  As(ja»)]  >  o  [T'Ujw)  -  I]  Vae[0,l]  (3-21) 

which  in  turn  will  be  satisfied  if  (using  (3-19)): 

o  [I  +  a  As(ju>)]  >  —  Vote[0,l]  (3-22) 

~  1-e 

After  multiplying  the  bracketted  term  on  the  left  by  (the  unitary  matrix)  # 
and  on  the  right  by  U,  condition  (3-22)  becomes: 

o  (I  +  a  Z  VHU]  >  ~  Vac [0,1]  (3-23) 

1-e 


Since  U  is  a  unitary  matrix  and  vi  has  unit  magnitude,  the  maximum  singular 
value  on  the  right  of  (3-25)  satisfies: 


Substituting  (3-26)  into  (3-25)  and  maximizing  the  right  side  over 


ue[0,l],  we  ns,  and  As  e  Vs  implies  that  (3-25)  will  be  satisfied  if 


1  +  a  0^  cos  $  e3®  0 


>  -  +  sup  [a^  sin  <j>]  =  r  Vae[0,l] 

e-1  weOc 


(3-27) 


Condition  (3-27)  will  be  satisfied  if 


1  +  a  oj  cos  <j>  eJ®|  )  r 


Vote  [0,1] 


(3-28) 


Assume  there  exists  and  ag  e  [0,1]  such  that  (3-28)  does  not  hold.  Then, 


since  Vs  is  convex  and  0  e  Vs, 


Aso  =  ag  As 


(3-29) 


is  an  element  of  Do  with 


As0  “1  *  °10  cos  4>0  e^G° 


c\  cos  <j>  eJ® 


(3-30) 


Hence 


1  +  ojg  cos  <|>o  e^ ®0  |  <  r 


which  contradicts  (3-18).  Thus  condition  (3-28)  must  hold.  Since  (3-28)  is 


sufficient  for  (3-5)  the  proof  is  complete. 


Theorem  3.4  provides  a  simple  graphical  test  for  the  robustness  of 


systems  with  rank  one  error  sources.  The  test  is  a  Nyquist-like  test  on  the 


(non-analytic)  function 


f(s)  =  ai  cos  <j>  eJ® 


u^  (s)  As  (s)  uj  (s)  Ag  e  Ds 


(3-31) 


\;a  .X’lr, 


where  ui(s)  is  the  first  column  of  U.  Condition  (3-18)  requires  the  Nyquist 
locus  of  each  f(s)  to  avoid  a  circle  in  the  complex  plane  of  radius  r  and 
centered  at  the  point  (-1,0)  (see  Fig.  3-2).  Specifically,  the  circle  Fr  is 
defined  as 


Fr  =  (x  e  <E:  |  x  -  (-1,0)  |  <  r}  (3-32) 

Since  Ps  is  a  convex  set,  any  point  inside  a  locus  for  a  particular  A8  e  V8 
(i.e.,  any  point  on  a  line  segment  from  the  origin  to  the  locus)  will  lie 
on  a  locus  for  some  other  As  e  Va.  This  means  that  the  set  of  all  points  on 
a  locus  for  some  As  e  Vs 

L  *  (x  e  I:  x  -  f(j<u)  for  some  As  e  Vs,  u>  e  fls}  (3-33) 

is  a  solid  region  in  the  complex  plane  (as  illustrated  by  the  shaded  region  in 
Fig.  3-2).  Hence,  the  region  can  be  completely  characterized  by  its  boundary. 

Thus,  the  test  for  stability  can  be  conducted  by  the  following  steps: 

1.  Determine  r  (from  (3-17)) 

2.  Draw  the  circle  Fr  (see  (3-32)) 

3.  Plot  the  boundary  of  L  (see  (3-33)) 

4.  Check  that  F r(1  l  ■ 

Theorem  3.4  is  reminiscent  of  the  various  circle  theorems  [19]-[21]  that 
are  available  for  nonlinear,  time  varying  and  multivariable  systems.  However, 
this  result  incorporates  the  error  source  phase  and  directionality  structure, 
as  well  as  the  usual  magnitude  variation  of  the  error  source.  The  ability  to 
incorporate  phase  structure  is  also  more  general  than  the  structures  that  can 
be  analyzed  using  the  structural  singular  value  [22 ]— [ 23 ] . 


Figure  3-2.  Graphical  Test  for  Rank  One  Error  Sources 


The  set  Vs  defines  the  structural  information  about  the  error  source  that 
is  required  by  Theorem  3.4.  While  the  rank  one  structure  is  imbedded  in  the 
theorem  statement,  VB  must  still  define  the  phase  and  directionality  informa¬ 
tion.  This  information  can  be  specified  either  by  the  complex  boundary  func¬ 
tion  f(joi)  and  the  second  term  of  (3-17),  or  by  the  real  functions  ai(ju)), 
<j»(ju>)  and  Q(ja>).  < 

Directionality  information  is  incorporated  through  the  second  term  on  the 
right  side  of  (3-17).  This  terra  is  a  bound  on  the  how  much  the  loop  coupling 
can  be  modified  in  the  nominal  system.  If  the  relative  change  between  the 
nominal  and  true  system  input-output  directions  is  large  in  a  frequency  range 
that  also  has  a  significant  magnitude  variation  then  this  term  will  be  large. 
This  bound  can  be  very  conservative  (a  result  of  (3-27)).  However,  the 


,1 

1 

1 


conservativeness  can  be  removed,  if  necessary,  at  the  expense  of  a  more 
complex  test  and  a  more  complex  characterization  of  the  set  t)s. 

Finally, -it  should  be  noted  that  Theorem  3.4  can  be  used  to  determine  how 
tight  the  approximation  (3-16)  must  be  for  a  given  error  source  set  Vs,  and 
how  much  stability  margin  is  guaranteed  for  the  true  system.  Given  a  set  t?s, 
define 


sup 

meQ 


|  oi  sin  <J>| 


s 


^8?  ^8 


(3-34) 


and 


*2 


inf 

men  e 


|  1  +  oi  cos  $  eJ®  | 


•s 

Ase^s 


(3-35) 


Then  if 


*1  +  42  <  1 


(3-36) 


there  is  a  value  of  e  such  that  Theorem  3.4  holds.  In  particular,  any  value  e 
that  satisfies 


c  < 


l-4l-A2 


2-4i-42 


(3-37) 


will  work.  This  value  in  turn  determines  how  tight  the  control  loop  must  be 
over  the  frequency  range  ft8,  The  controller  can  be  chosen  to  satisfy  (3-17) 
and  thus  guarantee  stability  (assuming,  of  course,  that  such  a  controller 
choice  does  not  violate  any  other  robustness  or  performance  constraints). 
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Once  a  value  for  e  has  been  determined,  the  graphical  test  can  be  used  to 
determine  the  guaranteed  stability  margin  m  of  the  closed  loop  system  with  the 


rank  one  perturbation  Included.  In  this  context,  the  guaranteed  stability 


margin  is  magnitude  of  the  largest  allowable  unstructured  error  in  either  the 


inverse  of  the  nominal  closed  loop  transfer  function  T”l(s)  or  in  the  error 

0 


source  characterization  Vs  that  can  be  tolerated.  The  value  m  is  given  by: 


m  =  inf  |  x-y  | 
xeFr 
yeL 


(3-38) 


That  is,  m  is  the  minimum  distance  from  Fr  to  L  (see  Fig.  3-2). 


3.3.4  Weak  Coupling  Errors 


The  third  principal  error  source  is  the  interconnection  dynamics  that 


result  from  connecting  the  subassemblies.  This  error  source  will  be 


characterized  by  the  set  of  interconnection  matrices  flw.  Let  {Pii0}N  denote 

i=l 


the  nominal  subsystem  models  of  the  decoupled  plant.  Then  the  nominal  model 


P0  =  diag  [PiioJ 


(3-39) 


and  the  error  between  Pq  and  any  true  model  P  is: 


p~l  p 
11  12 


•  •  P 

11  IN 


p“l  p 
22  21 


p“*  P 
22  2N 


P-1  p 


NN  N1 


p-1  p 


NN  N2 


.  0 


(3-40) 


Hence,  the  interconnection  error  set  is: 


=  [A  =  [Aij]:  A^£  =  0,  i=l,...,N} 


(3-41) 


A  number  of  systems  analyses  have  been  developed  to  handle  the  intercon¬ 


nected  systems  problem  (see  [24]  for  a  detailed  survey).  Most  have  assumed 


that  the  set  Vw  has  some  form  of  weak  interconnecting  structure.  The  most 


relevant  results  are  those  of  Bennett  and  Baras  [1],  and  Limebeer  [2].  Both 


utilize  the  concept  of  block  diagonal  dominance  [25]  to  prove  stability  of  an 


interconnected  system  given  a  nominally  stable  block  diagonal  system.  The 


following  Theorem  presents  a  similar  result  cast  in  the  framework  of  Theorems 


3. 1-3.2. 


Theorem  3.5:  (Weak  Interconnection  Error  Robustness  Test).  Assume  that 


assumptions  (A1)-(A6)  hold.  Let  (T0i£}N  denote  the  closed  loop  transfer 

i=l 


function  of  the  nominal  uncoupled  subsystem.  If  for  each  A  e  *>w 


0  f^oiil  ^ 


Li  • <4ij>] 


V  i=l , ...  ,N  V  a>  e  fly, 


(3-42) 


then  condition  J.5  holds  for  each  w  e  fy,. 


Proof:  Condition  (3.42)  can  be  rewritten  as: 


a  [T_l  ]  >  l  a  (Aij)  V  i=l,...,N  V  u  e  ^ 
oii  j-i 


(3-43) 


inequality  (3-43)  implies  that  (T~^  +  A)  is  block  diagonally  dominant  (see 

0 


[25]),  which  in  turn  implies  that  (T“^  +  A)  is  nonsingular  for  each  A  e  5W 


[25].  Thus 


o  (T“l  +  A)  >  0  V  A  e  Vw  Vu  e  !L 
”  0 


(3-44) 


Let  a  e  [0,1].  Then,  for  any  A  e  Vw,  aA  e  t?w.  By  (3-44), 


s 


a  (T-l  +  aA)  >  0  V  to  e  O, 

0  (3-45) 

Since  a  e  [0,1]  is  arbitrary,  (3-^5)  is  identical  to  (3-5). 

■ 

3.4  SUBSYSTEM  ERRORS  FOR  THE  DECOUPLED  SYSTEM 

The  first  error  source  results  from  the  unmodeled  higher  frequency  dynam¬ 
ics  due  to  the  reduced  order  design  model.  Since  such  neglected  dynamics  are 
almost  invariably  harmful,  it  is  appropriate  to  treat  them  as  unstructured 
multiplicative  perturbations  to  the  plant  as  represented  by  (3-1)  with  A(s)  e  Vm 
and  t>m  defined  as  in  (3-10).  Theorem  3.3  and  condition  (3-11)  give  the 
appropriate  robustness  test  for  errors  in  the  class  Vm, 

The  multiplicative  error  between  the  5  element  truth  model  (Pt(s))  and 
the  3  element  design  model  (Po(s))  defined  in  Section  2  is  given  by  (3-2).  A 
plot  of  the  inverse  of  the  largest  and  smallest  singular  values  is  shown  in 
Fig.  3-3.  The  bound  function  £m(ju))  can  be  taken  to  be 

4 

grrjd)  VOTTT  (3-46) 

Note  that  if  this  were  the  sole  error  source  a  bandwidth  of  approximately 
1-2  rad/ sec  could  be  achieved. 

3.5  SUBSYSTEM  ERRORS  FOR  THE  CONNECTED  SYSTEM 

The  second  error  source  results  from  the  dynamic  modification  of  the 
subsystem  transfer  function  caused  by  the  connection  of  the  two  subsystems. 

Let  Pt(s)  denote  the  subsystem  transfer  function  obtained  from  the  6  element 
analysis  model  of  the  coupled  system,  and  let  P()(s)  denote  the  3  element 
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1.76s  +  1 


Figure  3-3.  Inverse  Multiplicative  Error  Between  the  Nominal  3 
Element  Design  Model  and  the  5  Element  Truth  Model 

design  model  of  the  subsystem.  The  inverse  of  the  largest  singular  value  of 
the  relative  error  (3-2)  is  shown  in  Fig.  3—4.  Note  that  the  most  significant 
modification  to  the  model  is  for  frequencies  less  than  .1  radians/ second.  In 
fact,  if  this  error  were  completely  unstructured  condition  (3-2)  would  require 
the  closed  loop  transfer  function  to  have  a  magnitude  less  than  .2,  and  the 
control  objectives  would  not  be  attainable  within  this  control  architecture. 

Fortunately,  the  error  depicted  in  Fig.  3-4  is  highly  structured.  A  more 
detailed  and  accurate  analysis  can  be  performed  by  applying  the  rank  one 
analysis  tool  developed  in  Subsection  3.3.3. 

The  interconnection  of  two  subsystems  by  a  single  beam  can  be  viewed 
as  a  single-input  single-output  loop  closure.  Figure  3-5  depicts  a  single 
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MAGNITUDE 


FREQUENCY  (rad/sec)  r-2548 


Figure  3-4.  Inverse  Multiplicative  Error  Between  the  3  Element 

Subsystem  Design  Model  and  the  6  Element  Analysis  Model 
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Figure  3-5.  Subsystem  Modification  Structure 
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VvU 


subsystem  described  by: 


y 

yi 


Po 

poi 

pio 

pii 

u 

UI 

_  __ 

(3-47) 


where  u  is  the  2  dimensional  subsystem  input  vector,  y  is  the  2  dimensional 
subsystem  output  vector,  uj  is  the  1  dimensional  interconnection  input,  and  yj 
is  the  1  dimensional  interconnection  output.  The  second  subsystem  is  assumed 
to  have  the  interconnection  transfer  function 


ui  -  P  yi 


(3-48) 


Following  the  interconnection,  the  subsystem  transfer  function  from  u  to 
y  becomes: 


y  -  Ipq  +  Pox  -  p  Pnr1  ?  Piolu 


A  _ 

*  u 


(3-49) 


The  relative  error  between  Pt  and  Pq  is  given  by  (3-2) 


&s  “  P’1  tPoi  U  "  P  Pll)"1  P  Piol 


(3-50) 


Since  the  bracketted  term  on  the  right  side  of  (3-50)  has  rank  one  for  all 
values  of  s  (Poi  and  Pio  are  2x1  and  1x2  matrices,  respectively)  the 
perturbation  has  the  rank  one  structure  examined  in  subsection  3.3.3. 

The  set  Ve  of  all  subsystem  interconnection  errors  can  be  characterized 
by  a  boundary  function  f(jw)  whose  Nyquist  locus  is  the  boundary  of  the  set  of 
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all  Nyqulst  loci  of  &s  c  Ds«  For  the  purposes  of  this  report,  we  will  use  the 
relative  error  between  the  nominal  3  element  subsystem  design  model  and  the 
subsystem  transfer  function  generated  by  the  6  element  interconnected  analysis 
model  to  define  the  function  f(j<o).  Thus,  f(ju>)  will  be  defined  as: 

f(jw)  =  <*1  cos  <{>  eJ®  (3-51) 

where  oj,  <f>  and  0  are  defined  by  (3— 13 )— (3—1 5 )  applied  to  the  relative  error. 

The  Nyquist  plot  of  f(jus)  is  shown  in  Fig.  3-6.  This  locus  will  avoid 
any  circle  centered  at  the  point  (-1,0)  with  a  radius  less  than  0.6.  Hence, 
if  r  (as  defined  by  (3-17))  satisfies 

r  <  .6  (3-52) 

then  the  rank  one  robustness  test  (Theorem  3-4)  will  be  satisfied.  For  this 
relative  error 

sup  [oi  sin  $]  <  10“4 

u>eR  (3-53) 

and  its  contribution  to  r  is  negligible.  As  long  as  the  control  system  design 
satisfies  (3-16)  with 

,  e  <  .38 

over  the  frequency  range  for  which  the  error  source  is  dominant 

08  =  { u)  <  .2}  (3-54) 

the  test  will  still  be  satisfied. 
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REAL  R-2M9 

Figure  3-6.  Nyquist  Plot  of  oj  cos  <j>  eJ® 

Additional  information  regarding  the  behavior  of  the  coupled  decentral¬ 
ized  system  can  be  readily  obtained  from  (3.39)  once  a  value  for  e  has  been 
specified.  For  e=0,  the  guaranteed  margin  is  plotted  in  Fig.  3-7  versus 
frequency.  Based  on  this  margin,  we  would  expect  the  coupled  system  using 
decentralized  control  system  designs  to  have  a  mode  with  a  natural  frequency 
of  .1  rad/ sec  and  a  peak  amplitude  of  approximately  6  db.  Since  the  actual 
design  will  require  a  nonzero  value  for  e,  the  peak  will  be  somewhat  higher. 
For  e  *  .25  (i.e.,  Tq(s)  within  -2,5  db  of  the  identity  over  the  frequency 
range  08) ,  the  peak  would  be  9db. 

3.6  INTERCONNECTION  ERRORS 

The  final  error  source  is  the  interconnection  transfer  function  between 
the  subsystems.  This  error  class,  Qw,  is  characterized  by  the  largest 
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SUBSYSTEM  COUPLING  ERROR 


Figure  3-7.  Stability  Margin  of  Coupled  Analysis 
System  With  Structured  Error 


relative  norms  of  the  interconnection  transfer  functions  {fcwi(w)}N  .  The 

i*l 

class  K  for  the  two  subassembly  problem  of  Section  2  is  defined  as: 


'w 


(A(s) 


0 

_f2l(s) 


Ai2(s) 

0 


:  a  {Ay  (jaj)}  <  Awi(aj),  i-1,  2;  j*i} 

(3-55) 

The  relative  interconnection  error  will  be  determined  from  the  3  element 

design  model  (Po(s))  and  the  6  element  analysis  model  (Pt(s)).  A  plot  of  the 

inverse  error  a  [P"l  P  is  shown  in  Fig,  3-8  for  subsystem  1.  Due  to  the 
011  12 

symmetry  of  the  problem,  the  plot  of  the  error  a  [P“l  Poi for  subsystem 

022 

2  is  Identical,  As  can  be  seen  from  Fig.  3-8,  the  bounds  (^wi^i*l,2  can  be 
chosen  to  be: 


fcwi(u))  *  -20  dB  i*l,2,  iii  e  ^ 
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(3-56) 


Of  -  *  *; 


.  V  %  *. 


INVERSE  RELATIVE  COUPLINC  ERROR 


l 


! 

1 


Figure  3-8.  Inverse  of  the  Relative  Interconnection  Error 

where  flw  is  the  frequency  range: 

ft*  -  {w:  .5  <  a)  <  2}  (3-57) 

Theorem  3.5  Indicates  that  this  error  will  not  be  significant  as  long  as  the 
nominal  subsystem  closed  loop  transfer  functions  are  small  («  20db)  in  the 
frequency  range  ftw. 

3.7  SUMMARY  OF  THE  ERROR  ANALYSES 

I 

The  analyses  of  the  error  sources  in  subsections  3. 4-3. 6  can  be  combined 
to  provide  specifications  for  the  decentralized  control  system  designs.  The 
total  error  A  can  be  written  in  terms  of  the  three  principal  errors: 

A(s)  -  Aa(s)  +  As(s)  +  Aw(s)  (3-58) 
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Each  of  the  three  error  sources  is  dominant  in  a  different  frequency  region. 

The  latter  two  (As  and  Aw)  possess  significant  structure  (rank  one  error  and 
weak  interconnection,  respectively)  over  the  frequency  region  in  which  they  v 
are  dominant.  Each  of  these  error  sources  (As  and  Aw)  is  relatively  small 
outside  the  frequency  region.  They  can  be  regarded  as  unstructured  errors 
outside  the  dominant  regions,  and  can  be  combined  with  the  error  source  Am 
with  insignificant  effects. 

Thus  we  can  decompose  our  stability  analysis  and  design  specification 
into  three  frequency  regions: 

=  {oj:  a)  <  .2} 

Sly  *  (ai:  .5  <  a)  <  2} 

«m  -  * 

The  design  specifications  can  then  be  derived  from  the  analyses  of 
subsections  3. 4-3. 6  and  Theorems  3. 3-3.5.  Theorem  3.4  implies  that  the 
transfer  function  To(s)  for  each  subsystem  must  be  nearly  unity  in  the 
frequency  range  Jl8.  That  is, 

To(jw)  «  1  id  e  «s  (3-59) 

The  weak  coupling  Theorem  (Theorem  3.5)  implies  that  nominal  transfer  function 
design  must  be  less  than  20dB  to  avoid  the  interconnection  error  in  the 
frequency  range  flw.  There  is  actually  a  more  severe  constraint  in  this  region 
resulting  from  the  multiplicative  error  source  Ara.  This  source  consists  of 
the  nominal  subsystem  modeling  error  (Fig.  3-3),  and  the  subsystem  connection 
error  (Fig.  3-4).  If  these  are  assumed  to  combine  in  the  most  harmful  manner, 


the  total  inverse  error  would  be  bounded  by 


o  (A(jaO)-l  >1.9  u  e  !i„  (3-60) 

Hence,  by  Theorem  3.5  the  nominal  subsystem  transfer  function  designs  must 
satisfy 

7  [To(jw)]  <1.9  (3-61) 

Finally,  the  unstructured  multiplicative  error  source  becomes  significant 
at  frequencies  greater  than  2  radians/sec.  From  Fig.  3-3,  the  bound  Jtra( j uj) 
(Eq.  (3-46)),  and  Theorem  3-3,  the  nominal  subsystem  transfer  function  must 
satisfy. 

4 

0  tTo(jw)]  <  r  Vu> 

/3.1  0)2  +  1  (3-62) 

Conditions  (3.61),  (3.63)  and  (3.64)  provide  specifications  for  the 
subsystem  designs  that  will  guarantee  that  the  closed  loop  system  will  be 
stable,  both  while  decoupled  and  while  connected. 
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SECTION  4 


DECENTRALIZED  CONTROL  SYSTEM  DESIGN  AND  ANALYSIS 

4.1  LINEAR  QUADRATIC  REGULATOR  DESIGN 

A  control  system  was  designed  for  each  of  the  subassemblies  to  satisfy 
the  decentralized  design  requirements  developed  in  Section  3.  To  recapitu¬ 
late,  these  requirements  were: 


T0  (jm)  «  1 

w  <  .2 

(4-1) 

#  T0  (jui)  #  <  1.9 

1  <  a)  <  2 

(4-2) 

4 

II  Tq  ( jui)  11  <  - 

Vrn 

(4-3) 

/3.1m2  +1 

The  design  goal  is  to  obtain  the  largest  bandwidth  possible  (in  the  sense  de¬ 
scribed  in  subsection  3.1)  with  zero  steady  state  errors  in  the  mass 
separation. 

The  subassembly  control  system  designs  were  developed  using  the  Linear 
Quadratic-Gaussian/Loop  Transfer  Function  Recovery  (LQG/LTR)  procedure  [13]. 
Integral  states  were  appended  to  produce  the  desired  zero  steady  state  errors. 
The  subassembly  design  models  can  be  represented  in  state  space  form  as: 

*m  *  +  Bnu  (4-4) 

Define  the  integral  states: 

e  =  y  (4-6) 
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(4-7) 


(4-8) 


(4-9) 


(4-10) 

Define  the  regulator  cost 


where 

J  -  f  Qd  z  +  u^  Ru  dt 

0 

(4-11) 

Qd  *  cDT  cD 

(4-12) 

cd  -  [ci  c7|  * 

(4-13) 

,  R  -  Pi 

(4-14) 

The  the  LQ  regulator 

solution  is  given  by 

u  =  "  zm  =  “  Gp]  |~e  “I 

w 

(4-15) 

Gra  =  R_1  Bdt  K 

(4-16) 
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,  Then  the  regulator  design  system  becomes: 

i 

i 

! 

j  z  =  Ajjz  +  Bj)U 


where 


where  K  is  the  solution  of  the  regulator  Riccati  equation 


Adt  K  +  K  AD  +  QD  -  K  BD  R_1  Bdt  K  =  0 


(4-17) 


The  loop  transfer  function  L(s)  of  the  LQ  regulator  can  be  approximated  by 


L0(s)  =  [-  C!  Cm  (sI-Am)“l  Bm  +  Cp  (sI-Am)-l  Bm] 


(4-18) 


Equation  (4-18)  was  used  to  specify  Cj  and  Cp  to  satisfy  the  basic 
design  goals  summarized  at  the  beginning  of  this  section.  The  integral  weight 
Cj  was  chosen  to  satisfy 


1 


Lq(s)  -  -  I 
s 


s|«l 


(4-19) 


This  was  achieved  by  specifying 


Cl  =  {Cm  Am  *  Bm]  1 


(4-20) 


The  state  weighting  matrix  Cp  was  chosen  to  penalize  only  the  outputs: 


Cp  -  X  c* 


(4-12) 


The  remaining  degrees  of  freedom  (i.e.,  the  choice  of  X)  were  used  to  match 
the  high  frequency  singular  values  of  (4-18): 


X  ”  [Cm  Am  ^  Bm]  1 


(4-22) 


Finally,  the  scale  factor  p  of  the  input  penalty  matrix  was  chosen  to 
maximize  the  bandwidth  of  the  LQ  regulator  by  choosing  its  value  as  small  as 
possible  while  still  satisfying  (4-2).  The  resulting  value  of  p  was: 


p  =  .1 
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(4-23) 


The  resulting  LQ  regulator  design  has  closed  loop  eigenvalues: 


-.301  ±  1.53j 
-.223  ±  1.01  j 
-.280 
-.292 


subsystem  modes 


integrator  poles 


(4-24) 


The  singular  values  of  the  closed  loop  system  transfer  function  from  system 
input  to  controller  output  is  shown  in  Fig.  4-1.  The  singular  values  of  the 
return  difference  at  the  system  input  is  shown  in  Fig.  4-2. 


FREQUENCY  (r«d/sec) 


*-2552 


Figure  4-1.  Linear  Quadratic  Regulator  Transfer  Function 


The  primary  features  of  note  with  respect  to  the  LQ  designs  are  that 
the  design  constraints  (4-l)-(4-3)  imposed  by  the  error  analysis  of  Section  3 
are  satisfied.  However,  both  condition  (4-1)  and  (4-3)  are  barely  satisfied 
by  this  design.  Consequently,  we  will  expect  that  the  unraodeled  mode  that 

Vs 

determined  (4-3)  is  likely  to  be  very  oscillatory  when  this  controller  is 
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0.01  0.1  1.  10. 

FREQUENCY  (rad/sec)  R-2S53 


Figure  4-2.  Singular  Values  of  the  Linear  Quadratic  Regulator  Return 
Difference 

applied  to  the  truth  system.  Also,  the  coupling  pole  that  is  expected  to 
occur  at  s  *  ±  .lj  will  be  somewhat  more  lightly  damped  than  the  analysis 
based  on  the  approximation  (4-1)  would  indicate.  In  fact,  the  amount  by  which 
(4-1)  is  violated  results  in  a  value  for  e  (see  (3-14))  of  .25.  As  subsection 
3-5  indicates,  we  should  expect  an  approximately  9db  peak  in  the  closed  loop 
transfer  function  of  the  true  system. 

Finally,  we  should  emphasize  that  this  design  was  developed  specifically 
to  test  the  a  priori  analysis  procedures  used  in  Section  3.  There  are^  numer¬ 
ous  features  of  the  preceding  design  that  could  prove  undesirable  for  a  true 
system.  However,  this  design  will  serve  the  purpose  of  pushing  the  design  to 
the  limits  of  design  specifications. 


4.2  LOOP  TRANSFER  FUNCTION  RECOVERY 

The  loop  transfer  function  recovery  procedure  developed  by  Stein  and 
Doyle  [12]  and  Kwakernaak  [26]  was  applied  to  the  model  system  (4— 4) — (4—5)  to 
produce  an  output  feedback  compensator.  The  Kalman-Bucy  filter  design  used 
the  spectral  density  matrices 

-  B„  Bj  (4-25) 

0  =  51  (4-26) 

to  produce  a  filter  design: 

xm  =  4jn  xm  +  Bmu  +  Hm  (y  -  (^  xm)  (4-27) 

where 

^  »  ICj  Q”1  (4-28) 

and  l  is  the  solution  of  the  filter  Riccati  equation: 

A«E  +  SAj  +  2m  -E  Cj  0-1  CmE  =  0  (4-29) 

The  LQG  feedback  compensator  that  results  from  the  Kalman-Bucy  filter  (4-27) 
and  the  LQ  regulator  design  (4-7),  (4-15)  is  described  by  the  equations: 

xm  =  -  Bm  Gp Hm  Cm)  xm  -  Bm  Gje  +  Hmy  (4-30) 

e  *  y  (4-31) 

u  =  “Gp  ~  Gje  (4-32) 

As  6-*-0,  the  loop  transfer  function  of  the  overall  closed  loop  system 
(evaluated  at  the  plan  input)  approaches  the  loop  transfer  function  of  the  LQ 
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regulator.  For  the  LQ  design  of  subsection  4.1  and  the  design  constraints 
(4-l)-(4-3),  the  choice 

6  =  10-5  (4-33) 

resulted  In  a  satisfactory  approximation  to  the  LQ  loop  transfer  function 
(accurate  for  frequencies  less  than  8  radians/second).  The  KBF  eigenvalues 
for  this  value  of  6  are 


-9.25  ±  9.36j 

-5.73  ±  5. 83 j  (4-34) 

The  singular  values  of  the  closed  loop  transfer  function  for  the  design  model 
(4-4)-(4-5)  and  controller  defined  by  (4— 30)— (4— 32 )  are  shown  in  Fig.  4-3. 

The  same  controller  was  applied  to  the  subsystem  truth  model.  The  singular 
values  of  the  resulting  closed  loop  transfer  functions  are  shown  in  Fig.  4-4. 
As  expected,  the  mode  corresponding  to  the  unraodeled  dynamics  is  highly 
oscillatory  and  produces  a  small  gain  margin.  The  eigenvalues  of  the  closed 
loop  system  for  the  5  element  truth  model  are: 


-9.38  ±  9.48j  ) 
-5.81  ±  5.88j  ) 

-.03  ±  4.89j 
-.048  ±  4,65j  j 
-.236  ±  l,55j  \ 

-.172  ±  1.01 j  j 
-0.209  ±  ,0053j  | 


filter  eigenvalues 


unmodeled  dynamics 


subsystem  modes 
integrators 
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MAGNITUDE  (db) 


FREQUENCY  (rad/ sec) 


R-25S4 


Figure 


4-3.  Singular  Values  of  the  Closed  Loop  LQG  Design  Subsystem 


FREQUENCY  (rad/sec)  r-2555 


Figure  4-4.  Singular  Values  of  the  Closed  Lpop  LQG  Truth  Subsystem 


4.3  ANALYSIS  OF  THE  INTERCONNECTED  DECENTRALIZED  SYSTEM 


The  decentralized  subsystem  compensators  (one  for  each  subsystem)  de¬ 
veloped  in  subsection  4.2  were  applied  to  the  20th  order  truth  model  of  the 
coupled  system.  The  singular  values  of  the  transfer  function  of  the  closed 
loop  system  are  shown  in  Fig.  4-5,  and  the  closed  loop  eigenvalues  of  the 
resulting  32nd  order  system  are: 


-9.38  ±  9.48j  \ 
-5.80  ±  5.88j  ) 
-9.38  ±  9.48j  \ 
-5.80  ±.  5.88j  J 
-.03  ±  4 . 76 j  \ 
-.05  ±  4.54j  ) 
-.03  ±  4.7 6 j  \ 
-.05  ±  4.54j  ) 
-.03  ±  4.26  j 


subsystem  1  filter 

subsystem  2  filter 

subsystem  1  unmodeled 
dynamics 

subsystem  2  unmodeled 
dynamics 

coupling  dynamics 
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-.226  ±  1.46j  ) 
-.167  ±  . 953 j  ) 
-.225  ±  1.46J  ) 
-.166  ±  • 953 j  J 
-.0083  ±  ,115j  j 
-.272  ) 
-.228  ) 
-.272  \ 
-.228  ) 


subsystem  1  modes 


subsystem  2  modes 
coupling  mode 

subsystem  1 
integrators 

subsystem  2 
integrators 


Figure  4-5.  Singular  Values  of  the  Closed  Loop  Interconnected  Truth 
System 

The  qualitative  and  quantitative  features  of  Fig.  4-5  and  the  closed  loop 
eigenvalues  (4-36)  illustrate  the  validity  of  a  priori  analysis  of  Section  3 
and  the  resulting  design  specifications  (4—1 )— (4— 3) .  The  closed  loop  system 
is  stable,  although  only  marginally  as  a  consequence  of  the  deliberate  attempt 
to  push  the  design  specifications.  The  higher  frequency  unmodeled  dynamics 
produce  a  small  gain  margin  at  approximately  4  radians/ second.  The  low 
frequency  oscillation  between  the  subsystems  has  a  natural  frequency  of  .11 
radians/second  as  compared  to  the  predicted  value  of  .1  radians/second  from 
subsection  3.2.  Also,  the  peak  magnitude  of  the  transfer  function  due  to  this 
mode  is  about  9db,  which  is  precisely  the  prediction  of  subsection  3.5. 
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4.4  SECOND  LEVEL  CONTFOL  PROBLEM 

The  asymptotic  analysis  of  Appendices  A  and  B  indicate  that  the  second 
level  control  problem  can  be  approximated  by  a  system  of  two  masses  (each  with 
an  equivalent  value  of  the  total  mass  of  one  subsystem)  connected  by  a  beam 
with  a  longitudinal  wave  velocity  of  C2>  For  the  system  described  in  Section 
2  (Eq.  2-16),  a  2  element  finite  element  approximation  yields  a  system  with  an 
oscillatory  mode  of  approximately  .15  radians/second.  This  figure  agrees 
nicely  with  the  observed  mode  from  subsection  4.4.  However,  an  analysis  of 
the  zeroes  of  the  truth  system  and  decentralized  controller  from  the  second 
level  input  to  the  second  level  output  indicates  a  non  minimum  phase  zero  pair  at 

.049  ±  .63 j 

This  zero  pair  is  a  manifestation  of  the  fact  that  the  low  frequency  approxi¬ 
mation  that  results  from  Appendices  A  and  B  are  only  valid  for  a  limited 
frequency  range.  The  pair  limits  the  ability  to  extend  the  bandwidth  of  the 
closed  loop  system  without  more  extensive  modeling  of  the  overall  system. 

In  recognition  of  this  limitation  a  simple  low  gain  controller  was  de¬ 
signed  to  slightly  increase  the  damping  ratio  of  the  lightly  damped  low  fre¬ 
quency  coupling  mode.  The  controller  consists  of  a  simple  lead-lag  filter 
with  a  lead  break  at  .1  radians/ second  and  a  lag  break  at  1  radian/ second. 

The  magnitude  of  the  return  difference  of  the  second  level  closed  loop  system 

I 

is  shown  in  Fig.  4-6,  and  indicates  a  30  db  reduction  in  sensitivity  at  .11 
radians/ second.  The  magnitude  of  the  transfer  function  of  the  second  level 
closed  loop  system  is  shown  in  Fig.  4-7.  This  figure  Illustrates  the  reduc¬ 
tion  in  the  peak  of  the  coupling  mode.  The  closed  loop  eigenvalues  corres¬ 
ponding  to  this  coupling  mode  are: 

-.063  ±  .318 j 
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SECTION  5 


SUMMARY 

This  report  has  initiated  a  study  of  the  problem  of  robust  decentralized 
control  system  design  for  large  space  structures  that  will  be  assembled 
sequentially.  The  salient  features  of  this  problem  are  the  general  system 

v  ' 

characteristics  of  a  LSS,  the  requirement  that  the  control  system  be  com¬ 
pletely  decentralized,  and  the  requirement  that  the  decentralized  control  sys- 

.  *  ? 

tern  must  stabilize  the  total  structure  as  well  as  a  partially  assembled 
structure. 

All  of  the  relatively  few  studies  of  robust  decentralized  control  have 
concentrated  on  theoretical  studies  of  generic  decentralized  systems,  with 
only  mild  assumptions  about  the  structure  of  the  underlying  system.  While 
such  results  are  important  (and  have  been  used  in  this  report),  they  can  be 
conservative  when  a  problem  possesses  more  specific  system  structures.  The 
research  presented  in  this  report  has  identified  and  used  system  and  error 
structures  that  will  be  generically  present  in  the  problem  of  decentralized 
control  system  design  for  sequentially  assembled  large  space  structure.  The 
system  structural  analysis  was  conducted  qualitatively  using  the  infinite 
dimensional  singular  perturbation  analysis  of  the  Appendices  and  [8]— [9] . 

The  system  structures  were  evaluated  for  their  ability  to  support  a  completely 
decentralized  control  architecture  and  one  such  structure  was  chosen  for 
detailed  analysis  and  design. 
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The  conceptual  methodology  for  robust  decentralized  design  that  was  pro¬ 
posed  In  this  report  requires  the  identification  of  all  significant  error 
sources  and  performance  measures,  the  quantification  of  the  possible  detri¬ 
mental  effects  of  the  errors  on  the  decentralized  design,  and  the  decentral¬ 
ized  synthesis  of  a  control  system  that  minimizes  the  detrimental  effects. 
These  dominant  error  sources  were  identified:  unstructured  modeling  error, 
subsystem  modeling  error  induced  by  the  interconnection  of  the  subsystems,  and 
the  neglected  Interconnection  dynamics.  Appropriate  analysis  tools  were 
identified  or  developed  for  each  error  source.  Both  the  unstructured  error 
and  neglected  Interconnection  dynamics  were  quantified  using  existing  results. 
The  subsystem  errors  Induced  by  the  interconnection  were  sufficiently  large 
that  the  structure  of  the  error  had  to  be  exploited.  This  structure  was 
identified  and  a  new  robustness  analysis  technique  was  developed  to  produce  a 
less  conservative  quantification.  The  three  separate  analyses  were  then 
combined  using  the  integration  theorem  that  was  derived  in  Section  3. 

These  analysis  procedures  were  applied  to  a  specific  example.  The 
analysis  results  were  used  to  specify  design  constraints  for  the  decentralized 
synthesis  problems.  A  decentralized  design  was  then  developed  with  the 
specific  objective  of  testing  the  ability  of  the  analysis  results  to  predict 
the  properties  of  the  true  closed  loop  system.  The  properties  predicted  by 
each  of  the  individual  analyses  and  the  integrated  analysis  were  accurate. 

Several  extensions  to  this  work  should  be  pursued.  First,  the  design 
methodology  and  analysis  tools  have  only  been  applied  to  a  single,  simplified 
problem.  A  more  complete  problem  which  incorporates  more  realistic  models, 
error  sources,  and  control  objectives  should  be  addressed.  Secondly,  only  the 
completely  decentralized  control  architecture  has  been  studied  in  this  report. 
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As  noted  in  Section  2,  this  architecture  requires  a  specific  system  structure. 
Other  control  architectures,  such  as  a  partially  overlapping  control  system, 
may  be  useful  for  other  system  structures  and  for  relatively  homogeneous 
systems.  The  issues  raised  by  the  interaction  of  such  system  structures  and 
control  architectures  need  to  be  examined.  Finally,  the  possibility  of  using 
hierarchical  designs  for  enhancing  the  control  system  performance  following 
the  subsystem  assembly  is  promising  and  warrants  further  study. 
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ABSTRACT 


This  paper  examines  the  spectral  decomposition  of  a  class  of  operators 
that  depend  on  a  small  parameter  e.  The  convergence  of  the  eigenvalue- 
eigenvector  pairs  as  e+0  of  these  "stiff”  operators,  using  the  terminology  of 
[7]*,  is  investigated  with  the  objective  of  sheddding  light  on  the  singular 
behavior. 


^References  are  indicated  by  numbers  in  square  brackets,  the  list  appears  at 
the  end  of  this  Appendix. 
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SECTION  I 


INTRODUCTION 

This  paper  considers  the  eigenvalue  problem  of  the  so-called  "stiff" 
operators  using  the  terminology  of  [7].  Their  dependence  on  a  small  parameter 
e  suggests  some  spectral  "separation".  This  last  intuitive  idea  is  fully  in¬ 
vestigated  in  the  sequel. 

Many  physical  problems  can  be  described  by  models  containing  stiff 
operators.  Examples  of  such  problems  in  distributed  parameter  systems  are 
numerous.  Without  being  exhuastive,  the  examples  include  the  following: 

1.  Interfaced  media  having  very  different  properties  from  one 
region  to  another,  such  as  diffusivity  in  heat  conduction  or 
permittivity  in  electromagnetic  wave  propagation  [3]. 

2.  Continuous  stochastic  problems  when  the  noise  intensity  level 
is  different  from  one  part  to  another  of  a  medium  [2], 

There  are  several  motivations  for  studying  the  present  eigenvalue  prob¬ 
lem.  First,  the  "singular"  behavior  of  stiff  operators  can  be  better  under¬ 
stood  by  investigating  thier  spectral  decomposition  as  a  function  of  e. 

Second,  it  is  well-known  that  the'  solution  of  a  boundary  value  problem  can  be 
expanded,  at  least  theoretically,  in  infinite  series  of  the  eigenvectors.  The 
presence  of  a  small  parameter  may  be  advantageous  because  the  aforementioned 
solution  may  also  be  expanded  in  fractional  powers  of  e,  depending  on  its  con¬ 
vergence  as  e+0.  Therefore,  it  is  important  to  analyze  the  convergence  as  e+0 
of  the  eigenvalue-eigenvector  pairs  of  stiff  operators. 


The  eigenvalue  problems  involving  several  perturbed  operators  has  been 
studied  in  the  literature  [4]-(6],  {8]  and  the  references  therein.  However, 
the  spectral  analysis  of  stiff  operators  seems  to  be  new.  In  the  sequel,  a 
general  formulation  of  the  eigenvalue  problem  of  stiff  operators  is  presented, 
using  bilinear  forms  to  avoid  cumbersome  and  complex  boundary  and  interface 
conditions.  The  properties  of  the  eigenvalue-eigenvector  pairs  as  functions 
of  e  are  analyzed. 

This  paper  is  organized  as  follows.  The  eigenvalue  problem  formulation 
of  a  class  of  stiff  operators  involving  two  bilinear  forms  is  presented  in 
Section  2.  In  Section  3,  the  convergence  of  the  eigenvalues  and  their 
corresponding  eigenvectors  as  e+0  is  investigated.  In  Section  4,  a  numerical 
analysis  of  the  operator  of  Example  3-2  is  undertaken  to  illustrate  some  of 
the  properties  discussed  in  Section  3.  Finally,  the  last  section  contains 
some  concluding  remarks. 


SECTION  2 


EIGENVALUE  PROBLEM  FORMULATION 

In  this  section,  the  eigenvalue  problem  of  class  of  stiff  operators  is 
formulated.  Let  V,  H  be  two  given  real  Hilbert  spaces  such  that  V  is  dense  in 
H  and 

Al)  the  injection  of  V  into  H  is  compact. 

Let  V*  denote  the  dual  space  of  V.  After  identifying  H  with  H*,  one  has 
A2)  VCHCV*. 

Let  ai(<j>,i|>),  i  *  0,1  be  two  forms  on  V  such  that  the  following  assumptions 
hold: 

A3)  a£(<j>,'l')  Is  bilinear,  symmetric  on  V 

A4)  ai(<J>,i|0  is  continuous  on  V,  i.e.,  there  exists  such  that 
aiU,*)  <  Si  «<MV  l**Vt  V^V,  V'l'CV 

A5)  ai(4,$)  >  aiPi($)^»  where  a£  >  0  and  pi(*)  is  a  continuous  semi-norm 

on  V 

A6)  po(4>)  +  Pi($)  is  a  norm  equivalent  to  l$ly 

A7)  a^ ( 4> , )  *  0  on  VjCV,  where  is  an  infinite-dimensional  subspace  of 
V,  i-0,1. 

A8)  If  t-*-L{)(<|>)  is  a  continuous  linear  form  on  V,  null  on  Vo,  there 
exists  $eV  (modulo  Vq)  such  that 


ao(t,<»)  -  Lq(i|0,  V$eV 


Lee  ae($,+),  a($,<|>)  be  defined  as 


t- 


| 


& 

g 

l 

.1 

i 


ae(<J>,<>)  =  ao(4>,iJ>)  +  eai(4>,<|>)  (2-1) 

a($,i}>)  =  ag(  +  a^  ( <|>,^)  (2-2) 

Now  some  important  remarks  clarifying  the  above  assumptions  and 
definitions  are  in  order: 

Remark  2.1: 

It  can  be  easily  seen  from  (A3-A6)  that  ae(<f>,<|t),  a( 4» , 4»)  are  bilinear, 
symmetric,  coercive  and  bounded  on  V.  In  particular,  for  sufficiently  small 
e,  they  cctiory 

2  2 

acidly  <  ae(<ji,<t>)  <  v#4>Bv  ,  V$eV  (2-3) 

2  2 
al<j»#y  <  a(<M)  <  vl<t>lv  ,  V<j>eV 

where  a  (resp.  v)  is  independent  of  c  and  depends  solely  on  oq,  (resp.  v>o» 
v^)  and  the  semi-norms  Pl(*),  i  *  0,1. 

Remark  2.2: 

The  bilinear  forms  ae($,^i),  a($,i|>)  define  self  ad  joint  operators  [1] 

Ad,  ACf(V;V*), 

i .  e « , 

aeU.t)  *  (Ae^'P)  .  V  ♦  e  V 
a($,<|<)  -  (A*,*),  V  e  V 

where  (*,  •)  denotes  the  duality  pairing  between  V  and  its  dual  V*.  From  the 
preceding  remarks,  one  concludes  that  the  spectra  of  Ae,  A  are  subsets  of  1R+, 
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consisting  only  of  the  point  spectrum  [1],  15],  [12].  The  eigenvalue  problem 

for  Ag  is,  then,  to  seek  (yk,  xkJ  6  x  V  such  that 

e  e 

AeXk  -  YkXk  •  (2-4) 

e  e  e 

The  equivalent  variational  formulation  is 

Mxk,$)  -  Yk( Xk »  ♦)  V^eV  (2-5) 

e  e  e 

Now,  some  well-known  facts  are  summarized  in: 

Proposition  2.1: 

If  (A1-A6)  hold,  then  there  exist  unique  sequences  (Yk)"  e  ]R+, 

e  k=l 

{xk}°°  e  V  such  that  (2-4)  (or,  equivalently,  (2-5))  is  satisfied, 
e  k=l 

Furthermore, 

1)  0  <  y1  <  Y2  <  lim  Yk  “  +00 

C  £  (^4-f-oo  £ 

2)  {xM“  is  a  complete  orthonormal  set  in  H,  i.e.,  in  particular 

e  k*l 

(xk»X^)  “  $k*  (Kronecker  delta) 
e  e 

3)  The  multiplicity  of  each  eigenvalue  is  finite 

4)  The  eigenvalues  satisfy  the  following  minlmax  formula: 

Yk  M  min  >  max  ae(x,x) 

E  WCV,  dim  W-k  X  e  V,  1X1  -1  (2-6) 

H 

Proof:  See  [1],  [5],  [12] 

Remark  2.3: 

It  is  noteworthy  to  mention  that  Prop.  2. 1  is  valid  for  any  positive 
value  of  e,  e.g.,  e-1.  In  this  case,  the  eigenvalue-eigenvector  pais.  of  A  is 
obtained. 
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SPECTRAL  ANALYSIS  OF  THE  OPERATOR  Ae 

This  section  starts  with  a  series  of  lemmas  which  characterize  the 
various  properties  of  the  spectrum  of  Ae.  Then  the  convergence  of  the 
eigenvalues  and  their  corresponding  eigenvectors  as  e-K)  is  stated  and  proved 
in  Theorem  3.1.  Some  typical  examples  are  given  at  the  end  of  this  section 
to  illustrate  the  ideas. advanced  in  the  course  of  the  present  analysis. 

One  way  to  gather  information  about  the  behavior  of  a  single  eigenvalue 
as  e+0,  is  to  bound  it  from  below  and  from  above  by  known  functions  of  e. 
This  task  is  accomplished  in: 


Lemma  3.1: 

For  sufficiently  small  positive  e,  the  following  estimate  holds: 


eyk  <  Yk  <  Yk 

e  (3-1) 

for  k  ■  1,2....,  where  (yM“  are  the  eigenvalues  of  the  operator  A,  i.e., 

k-1 

they  satisfy  > 

a()(Pk>$)  +  al(Pk»$)  *  Yk(pk»$)»  V  ♦  e  V 


Proof:  For  sufficiently  small  e,  one  has 


ea(M)  <  ae(*,*)  <  a(*,*)  ,  V  ♦  e  V  (3-2) 

Using  the  mlnlmax  characterization  of  eigenvalues  (2-6),  one  readily  deduces 
(3-1)  from  (3-2). 
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Now  an  upper  bound  for  the  eigenvector  norm  in  V  is  derived: 


Lemma  3.2: 

If  xk  is  any  normalized  eigenvector  of  Ae,  corresponding  to  yk,  i.e., 
e  e 

H  Xk  11 H  =  then 

e 

aelx^i^  <  yk  <  yk  (3-3) 

e  V  e 

for  k  *  1,2, . . .  . 


Proof:  For  each  k,  the  sequence  yk  is  bounded  by  Lemma  3.1. 

e 

(2.5)  to  get 


Let  <{>  =  xk  in 
e 


a  (xk,Xk)  *  Yk(xk.Xk) 

e  e  e  e  e  e 


*  yk  <  yk 

e 


Now  one  easily  gets  (3-3)  by  using  (2-3). 

At  this  point,  the  tools  necessary  for  finding  the  limits  of  the  eigen¬ 
values  are  available. 


Lemma  3.3: 

The  sequence  {yk}-  is  decomposable  into  two  subsequences  {yk}“ 

e  k-I  e  k-1' 

{yk}”  such  that,  for  each  k, 
e  k*l 

V 

lira  Xk  -  o  (3-4) 

e+0  e 

lira  yk  ■  yk  >  0  (3-5) 

e-K)  e  0 


k  >  1,2,...  . 


Proof:  The  following  three-step  contradiction  argument  is  used  to  ascertain 
the  above  lemma. 
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1)  Suppose  llm  yk  =  0,  k  =  1,2 .  Let  V  =  Vq0V^,  where  the  ortho- 

e+0  e  0 

gonallty  is  that  of  V.  Take  v  e  W  and  write  it  as 

0 


v  =  £  (v,x’c)xk 
k=l  e  e 


(3-6) 


with  this  choice  of  v,  the  following  inequality  holds: 


ae(v,v)  =  ao(v,v)  +  eai(v,v)  >  C 


(3-7) 


for  some  strictly  positive  constant  C,  which  is  independent  of  e. 

00 

Using  (3-6),  ae(v,v)  *  £  Yk(v,xk)^,  which  converges  to  zero  as  e-K), 

k=l  e  e 

contradicting  (3-7). 

2)  Suppose  lim  yk  *  yk  >  o,  k  »  1,2...  .  Select  v  e  Vq  and  write  it  as 
e+0  e  0 

in  (3-6).  It  is  clear  that  the  following  inequality  holds: 


ae(v,v)  -  £  yk(v,xk)2  >  C 

k=l  e  e 


(3-8) 


for  some  strictly  positive  constant  C,  which  is  independent  of  e.  However 

ae(v,v)  “  eai(v,v),  which  converges  to  zero  as  e-K),  contradicting  (3-8). 

3)  Suppose  lira  yk  ■  0,  for  k  *  1,2, ...,£,  without  loss  of  generality. 
e+0  e 


V£  -  span{xl,x2..-.*,X^}  • 

e  e  e 

Select  v  e  Vq  ^  and  repeat  ste)p  2,  to  conclude  that  Vjj,  is  infinit  dimension¬ 
al,  unless  Vq  is  degenerate  (i.e.,  finite-dimensional).  An  identical  argument 
can  be  advanced  to  contradict  the  possibility  that  (3-5)  is  true  for 
k*l,2,...,£  (i.  finite). 

Now  decompose  (yk,  xk)“  into  (Ak,$k}“  if  lim  yk  *  0  and  into 
e  e  k*l  e  e  k*l  e 

,  otherwise. 

E  E  k»l 


v.mvi  %  i\,vi 


Some  attention  must  be  focused  on  how  A*1  converges  to  zero,  k  =  1,2,, 

e 

Lemma  3.4 

The  sequence  A^  converges  to  zero  with  a  linear  rate  A*4,  i.e., 
c  1 

Ak  -  A ke  +  o(e)  ,  k  =  1,2,... 
e  l 

Proof:  By  Lemma  3.1,  one  may  assume,  without  loss  of  generality  ,  that 


Ak  =  xkev  +  0(e) 
e  1 


(3-9) 


with  v  e  (0,1).  In  order  to  complete  the  proof,  it  suffices  to  show  that  v=l. 
If  <J>k  is  a  normalized  eigenvector,  i.e.,  ($k,  $k)fl  -  1,  corresponding 


e  e 


to  Ak,  then 

e 


Ak  ”  ao($k,  $k)  +  eai($k,$k) 
e  e  e  e  e 


(3-10) 


from  which  one  observes  that 


a 0Uk,$k>  <  0(eV)  (3-11) 

e  e 


Therefore,  ao ( » 4»fc)  ♦  0  as  £  0,  which  implies 

£  £ 


p  ( $k)  +  0,  as  E  ♦  0 
0  e 


by  (A5).  Using  (A6),  one  has 


p  (*k)  >  C  (3-12) 

1  £ 


for  some  strictly  positive  constant  C,  which  is  independent  of  e. 


Suppose  that  v  <  1.  Then,  from  (3-9)-(3-12)  one  concludes  that 


A-ll 


•  w  -  *- 


».V  V-~»' 


■  ■■  V  »■- 


■■  ***u*>»'' 


(3-13) 


xk  =  X*cev  +  o(e) 
c  1 


=  ag^k^k)  +  o(e) 
e  e 


.  e 

From  this,  there  exists  an  element  of  V,  <JTK  -  -  such  that 

e  ev/2 


Xk  =  ao(rk,rk)  +  0(el-v)  . 

1  e  e 


However,  such  a  claim  is  false  because 


(♦-k»$-k)H  =  -  +  +»  as  e  +  0. 


e  e  v 

e 


Since  the  injection  of  V  into  H  is  continuous,  ll<j»-k|lv  +  +®  as  e  +  0.  In 

e 

conclusion,  there  is  no  element  of  V  such  that  (3-12)  is  satisfied. 

Remark  3.1: 

For  v  <  1,  the  major  contribution  to  Xk  is  supplied  from  V1  (c.f.  (3-13)), 

e  0 

but  the  norm  of  the  contribution  is  concentrated  on  Vq  (c.f.  (3-12)),  which  is 
the  paradox. 

Hereafter,  the  focus  will  be  on  the  asymptotic  behavior  of  the  eigen¬ 
vectors.  The  following  lemma  summarizes  the  norm  bounds  of  the  eigenvectors: 
Lemma  3.5: 

Let  {X^,^^}00  ,  {1)^,^}“  be  as  in  the  proof  of  Lemma  3.3,  with  the 

e  e  k«l  e  e  k*l 

eigenvectors  normalized  in  H.  Then,  for  sufficiently  small  e,  the  following 
estimates  hold: 


1)  l<t>k|v  <  Ci 
e 


(3-14) 
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0-15) 


2)  /£"  li|tk#v  <  c2 

e 

k  *  1,2,... 

where  Cj ,  C2  denote  constant  Independent  of  e. 

Proof:  Use  Lemmas  3. 2-3.4. 

The  forthcoming  theorem  is  the  main  result  of  this  section.  It  states 
the  convergence  of  the  eigenvalues  and  their  corresponding  eigenvectors  as 
e+0. 

Theorem  3. 1 

Let  {Yk}°°  be  the 'eigenvalues  of  Ae  and  (x^l*  the  corresponding 
e  k-1  e  k=l 

normalized  system  of  eigenvectors.  Then,  given  a  sequence  of  e  converging  to 

zero,  {Yk,xk}“  can  be  decomposed  into  two  subsequences  {Xk,$k}“  and 
e  e  k=l  e  e  k-1 

(pk,<f>k}°°  which  have  the  following  asymptotic  properties,  for  each  k: 
e  e  k=l 

1)  Xk  0  linearly  in  e,  $k  *  $k  weakly  in  V 

e  e 

2)  pk  +  yk  >  0,  i|)k  +  \J^  weakly  in  H 

e  0  e 

where  {<j>k}“  and  {ij>k}"  satisfy 

k-1  k-1 


ai(4»k,x)  -  Xk(*k,X),  *k  e  V0CV,  VX  e  V0 

A  V 

(3-16) 

ao(+k»x)  *  Mk(4»k>x),  tyk  e  HiCH,  Vx  e  V. 

0 

(3-17) 

Proof:  Using  the  fact  that  l$klH  -  1,  and  the  estimates  (3-11),  (3-14),  one 

e 

conludes,  that,  given  a  sequence  of  e  converging  to  zero,  $k  $k  weakly  in  V 
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(hence  strongly  in  H  by  compactness).  From  (3-11),  it  results  that  <J>k  e  Vq. 

By  Lemma  3.4,  Xk  is  asymptotically  equal  to  Xke.  Hence  (2-5)  degenerates  into 
e  1 

(3-16)  in  the  limit.  Now, 


Pk  =  ag(  »j/k,  ijik)  +  ea  (ij>k,i|»k) 
e  e  e  lee 


(3-18) 


for  normalized  to  1  in  H.  Since  pk  is  0(1),  using  (3-18)  and  the  minimax 
e  e 

characterization  of  eigenvalues,  i.e.,  (2.6),  it  results  that  aQ(i|>k,i|)k)  is 

.  e  e 

0(1)  and  ai(i|(k,i|;k)  is  0(  -)  or  equivalently 
e  e  e 


P  (*k)  -  0(1) 
0  e 


(3-19) 


1 

p  ((j)k)  =  0(— ). 

1  e  /“ 


(3-20) 


Hence,  the  estimate  (3-15)  is  right.  Therefore 


+  +eo  e  0. 


(3-21) 


Note  that  ^  also  satisfies 


aoC^.x)  +  Ea,(xk»X>  ■  ^(♦'Sxh  V  x  e  V 

c  !  r  c 


e  e 


(3-22) 


From  (3—19 ) — (3—20) ,  one  deduces  that  ao(\|>k,x)  is  bounded  as  e  ♦  0  and 
l  C 

al(vkiX)  is  0( - ).  Since  *  1*  taking  formally  the  limit  as  e  0  in 


(3-22)  yields: 


<T 


ao^.x)  *  Uk('J»k>x) ,  v  x  e  V 

0 


(3-23) 


where  tk  e  Hi  (a  subspace  of  H). 


,  1 
St 

•-V-\  w 


Now  consider  the  following  boundary  value  problem 

«H)(vk,x)  +  eai(vk,x)  =  Uk(tk,x)»  V  X  e  V 
e  c  0 

which  admits  a  unique  solution  vk  e  V  for  positive  values  of  e.  As  e  0, 

e 

vk  +  strongly  in  H. 
e 

Let  wk  =  i|>k  -  vk  ,  which  satisfies 
e  e  e 

aO^'X)  +  eaj(wk,x)  =  bk(i|»k  -  ^.x) 
e  e  e  e 

+  (uk  -  Mk)(4»k,x)»  V  x  e  V. 
e  0 

The  left-hand  side  of  this  equation  goes  to  zero  as  e  ♦  0,  implying 

(^k  -  ^>k(x)  +  0  as  e  +  0. 
e 

Therefore 

tpk  +  \j>k  weakly  in  H. 
e 

The  results  of  Theorem  3.1  hold  for  negative  powers  of  e  as  well.  To 
illustrate  this  claim,  let 

1 

M<J>»X)  ■  "  t»-lC4>*X>  +  boU.'l') 

G 

where  bi(<J>,x),  i***l,  0  satisfy  (A3-A7).  Then,  one  has: 

Corollary  3. 1 

Let  {vk}“  be  the  eigenvalues  and  CxM"  the  corresponding  normalized 
e  k=*L  e  k-1 

eigenvectors,  derived  from 

Mxk,*>  "  vk(  Xk>$) >  V  «J>  g  V 

G  G  G 
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then,  given  a  sequence  of  e  converging  to  zero,  {vk,xM”  can  be  decomposed 

£  e  k=l 

into  and  ,  with  the  following  properties: 

e  e  k=l  e  e  k=l 

1)  Xk  -►  \k  >  o,  <j>k  <|i  weakly  in  V 

e  0  e 

k 

U-l 

2)  pk  =  -  +  0(1),  xk  "*■  Xk  weakly  in  H 

e  e  e 

where  {<|>k}"  and  (xk}“  satisfy 
k=l  k=l 

bo  {<}>k,x)  =  \k  (<J>k>x),  <)>k  £  Vq,  V  x  e  V 
b-l(xk,X)  -  Pk  (<frk,X>,  4>k  e  Hlt  V  X  e  V 

Proof 

Let  ai(<|i,x)  =  ebj[_i(()),x) ,  1=0,1  and  $k  =  e  yk  in  Theorem  3.1  to  get  the 

e  e 

desired  results. 

Remark  3.2: 

The  weak  convergence  in  Theorem  3.1  cannot  be  improved  in  general.  This 
will  be  illustrated  by  Example  3.2. 

Remark  3.3: 

A  careful  examination  of  the  steps  of  the  analysis  undertaken  in  the 

present  section  yields  the  following  observation:  the  weak  limits  {<|>k}“  , 

k=l 

{i|»k}“  form  an  orthonormal  system  in  H.  This  remark  is  of  paramount  impor¬ 
tance  in  approximating  the  solution  of  boundary  value  problems  involving  the 
operator  Ag. 

Now  some  examples  are  given  as  concrete  illustrations  of  the  above 
abstract  results.  Only  operators  of  order  less  than  or  equal  to  four  are 
considered,  due  to  their  frequent  usage  in  modeling  of  physical  processes. 
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Let  ft  =  Qo  ^1  S  be  a  bounded  set  in  IRn  with  boundary  T  =  rgUri.  the 
manifold  S  denotes  the  interface  between  Aq  and  Aj ,  as  indicated  in  Fig.  3-1: 


Example  3.1:  A  second  order  operator 

Let  H  =  L2(A),  V  =  H1^) 

0 

n  3<j> 

=  E  / - dx,  i  =  0,1 

j=l  Ai  3xj  3xj 

then  (2. A)  becomes: 


&  Xk  -  Yk  Xk  on  Aq 
eO  e  eO 


eA  Xk  =  Yk  Xk  °n  ^1 
el  e  el 


Xk  “0,  Xk, 

e0  r0  el  -  0 


Xk  =  Xk 
eO  el 


(3-24) 


where:  A  stands  for  the  Laplacian  in  3Rn. 

v  is  the  unit  normal  on  T  or  S,  outward  relative  to  Aq. 
In  this  example,  (3-16)  becomes 


^  *  0  on  An 
0 

A«j,k  ,  \k  *k  on  jj 
111  1 


$k  =  0 ,  ^  =*0 

1  rx 


(3-25) 
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Figure  3-1.  Examples  of  Interfaced  Sets 


which  is  a  Dirichlet  eigenvalue  problem  for  the  Laplacian  operator  in  Q^.  The 
subspace  Vq  of  V  is: 

Vo  *  {  X  e  V  :  Xq  =  0,  X^  e  (fy)} 

The  conclusions  of  Prop.  2.1  are  applicable  in  this  case.  Hence  {$*}"  in  a 

1  k“l 

complete  orthonormal  system  in 


Equation  (3-17)  becomes 


-  A  l'k  =  yk  0n  Jin 
10  0 


on  0^ 


1  r0 


0,  -  -  0 

3v  s 


(3-26) 


JV  L>  >.3iL' v  X 


which  is  an  eigenvalue  problem  with  mixed  boundary  conditions  for  the 

Laplacian  operator  in  89*  Again,  the  conclusions  of  Prop.  2.1  are  applicable 

in  this  instance,  provided  the  interface  S  is  sufficiently  smooth.  Therefore, 

{tyk}"  form  an  orthonormal  system  in  L^(J2q).  It  is  noteworthy  to  observe 
°  k=l 

that  i|/k  =  0  because  ipk  must  be  orthogonal  (in  (ft))  to  £=1,2,...  . 

I 

The  subspace  H]^  of  H  is  then 

3x0 

Hi  =  (x  €  H:  Xq  e  HUftojTo),  — -  =  0  on  S,  X^  =  0}  . 

o  V 

Remark  3.4 

In  Example  3.1,  one  can  consider 

n  n  3ij>  3ty 

ai( 4> =  Z  Z  /  a*  -  -  dx,  i  ■  0,1 

j*l  k=l  ft*  Jk  3xj  Sxfe 


which  a*  satisfies 

jk 


1)  a1  e  Ckftt) 
jk 


2) 


a* 

kj 


3) 


n  n  ,  n  « 

Z  Z  a1  £<Ck  >“1  1  °i  >  °» 

j-I  k=l  jk  J  k-I  k 


V  ;  e  IRn ,  C  *  0. 


The  discussion  therein  remains  unchanged. 

Example  3.2 

Let  fto  *  (-1.0),  -  (0,1)  S  -  {0},  r0  -  {-1},  Ti  -  {!}.  Then  the 


solutions  of  (3-25)-(3-26)  become  respectively 


=  Ot*)2 

1 

$k  =  0 
0 

<j>k  =  JT stn  kitx 


^  =  ((2k— 1 )  */2)2 


VJ 

$ 
Ti 


$k  =  cost (2k- 1 )  *x/2] 


^k  =  o 


The  exact 


eigenvalues  ^  must  satisfy  the  following  transcendental  equation: 


cos 


V^" 


sin 


+  fe  sin 


cos  -  “  o. 

e 


(3-27) 


n  ls  not  possible  to  solve  for  1*  ss  a  function  of  s.  However,  for  some 

£ 

c  p  _  „  c  >  1/ (A £+1 )2 »  l  an  integer,  (w1  »  )  given  below  is 

sequences  or  z,  e.g#,  e  1  /v1**  >  *  eJt  ei 

an  exact  eigenvalue-eigenvector  pair 

tul  «(;)2 

el  2 


’fr1  0  -  cos  -  x 
e£  2 


lj|l  I  *  cos 
ei 


.  I 


(3-28) 
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,%<)  v,r.-  f -y-MMta-b  ic3  w 5 h » £ w >  » ~  * * v~ wx vx vx vk yy  vk  w> w> 


Observe  the  oscilLatory  behavior  of  t|/^  as  l  *  “.  This  is  the  manifesta- 

ell 

tion  of  its  weak  convergence  in  (L2^))  to  zero.  Furthermore,  this  example 

shows  that  are  not  analytic  functions  of  e.  Consequently,  the  task 

e  k=l 

of  finding  these  eigenvectors  (in  order  to  solve  boundary  value  problems) 

is  impossible.  The  alternative  is  to  use  the  weak  limits  of  {y^}00  and 

e  k=l 

to  "approximate"  the  solutions  of  boundary  value  problems  involving 
k=l 

stiff  operators.  This  is  will  be  the  subject  of  a  subsequent  paper  (see  [9]). 

1  b  i. 

Similarly,  for  =  — ,  k-1,2,...,  {X^  ,  }  given  below  is  an  exact 

k2  ek  ek 

eigenvalue-eigenvector  pair 


Xk  =  (it)2 
ek 

*  /e^  sin  Sit  kitx 

ekP 


m  sin  kxx 

ckl 


(3-29) 


Note  that  in  (3-29),  depends  on  the  index  of  the  eigenvalue.  Hence, 

one  cannot  let  go  to  zero  and  observe  what  happens  to  this 

eigenvalue-eigenvector  pair«  However,  the  following  conclusions  can  be  drawn. 

First,  the  results  of  Theorem  2.1  show  that  ^  +  0  in  the  present  example. 

e0 

This  is  certainly  reflected  in  (3-29)  by  the  presence  of  the  factor  /c£ 

multiplying  <j>k  .  Second,  the  presence  of  /evT  in  the  argument  of  the  sin 
ek0 

function  in  (3-29)  indicates  that  "flattening”  occurs  in  ^  as  e  +  0. 

e 

Example  3.3:  A  fourth  order  operator 

Let  iJ  =  L2(n),  V  -  HZ(Q) 

0 

ai($,'|')  3  /  dx,  i  *  0,1 
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then  (2.4)  becomes: 


R 


&2yk  =  yk^k  on  jv 

eO  e  eO 


eA^k  =  ykyk  on 
el  e  el 


3xk 

eO 


X"  a  -  ■  0  on  r0 

eO  3v 


3x 


el 


el 


3\> 


0  on  Tj. 


Xk  =  Xk.  » 


3xk  3xk 

Xe0  el 


eO  el  3v 


3v 


3Axk  3&xk 

4xk  “  e  ^X^  »  - — —  “  e  ■— — 

eO  el  3v  3v 


Equation  (3-16)  becomes 


*k  »  o  on  0 
0 


a  Ak*k  on  {j. 
111 


3<|>k 


a  -  «  0  on  Ft 

l  3v 


3*k 

*k  a  -  *  0  on  S 

1  3v 


\ 


on  S 


(3-30) 


This  is  a  Dirichlet  eigenvalue  problem  for  the  biharraonic  operator  in  fy, 
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WWvU IVWiA,-,  .V-v-V-T,./-  ...  -  . 

v  ■  **  -  r.  «<.*  *.„  H  J  V\,  W".  y  ,  y  .  .  , 

"  ”  j  -  •-  -  *  ->  «\  s*.  xV  f  *  s« r 


V  iVLT  ■Wvf 


Equation  (3-17)  becomes 


A^tf^  =  ykd/k  on  Ho 
0  0  0 


^  =  0  on 


^,k  3  -  =  0  on  Tq 

0  3v 


Ai|/k  =  -  -  0  on  S 

0  3v 


(3-32) 


Identical  comments  to  those  of  Examples  3. 1-3.2  can  be  made  here,  provided  the 
function  spaces  are  changed  to  reflect  the  increase  in  the  operator  order  from 
two  to  four  (as  seen  in  (3-30)). 


Remark  3.5; 

Some  of  the  properties  discussed  in  this  section  are  shared  by  many  other 
stiff  operators,  which  do  not  fit  in  the  axiomatization  of  Section  2  (c.f. 
[9]). 
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SECTION  4 


NUMERICAL  RESULTS 

In  this  section.  Example  3.2  is  analyzed  numerically,  using  the  Finite 

Element  Method,  to  supplement  the  analysis  undertaken  in  the  preceding 

section.  The  set  Q  *  (-1,1)  is  divided  into  N  equal  intervals  of  length 
2 

h  =  The  roof  functions  {$*■)  ^“1  111],  [13]  are  selected  as  a  basis  for  the 
N  h  i-1 

finite  dimensional  approximation  of  H^(0).  The  finite  dimensional  approxi- 

0 

nation  of  the  operator  Ag  can  be  written  in  matrix  form  as 


where  the  entries  of  M*1  are 


Ah  .  (Mh)  -ljch 
e  e 


(M^)  ■  /  dx 

n  h  h 

i,j  -  1,2 . N-l. 


(4-1) 


In  the  forthcoming  example,  the  matrix  M^  can  be  written  explicitly, 


l»6i 


4  1 


Mh  - 


h 

6 


4 

1  4  1 


0 


l  4 
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For  all  computer  runs,  N  Is  selected  so  that  an  accurate  plot  Is  obtained  (N 

is  Indicated  under  each  plot).  The  subroutine  RSG  from  EISPACK  with  single 

precision  is  used  to  compute  the  eigenvalues  and  eigenvectors  of  A*1. 

e 


Example  4.1:  (Refer  to  Example  3.2) 

In  this  example,  the  entries  of  K*1  are 

e 


where 


d*1  d$j 

(Kh)  =  /  a(x)  — —  — —  dx 
s  i.j  n  dx  dx 


i,j  m  1,2,... ,N-1 


a(x) 


1  if  x  e(-l ,0) 
e  if  x  e(0,l) 


For  N  even,  K*1  can  be  written  as 
e 


2  -1 

-1  2 


Kh  -  i 

e  h 


-1 

0 


-1 

1+e 


0 

-e 


-e  -2e 


2e  -e 
-e  2e 


N 

■*•  -  +  1 
2 


(4-2) 


(4-3) 


The  eigenvalues  pi  are  tabulated  in  Table  4-1  for  e  -  0.1, 

E  E  E 


0.04,  0.01,  0.001.  The  corresponding  eigenvectors  are  plotted  in  Figures 
(4-la-4-ld) ,  (4-2a-4-2d),  (4-3a-4-3d). 
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Figure  4-3d.  <J»1  for  e  »  0.001,  N  *  150 

e 


TABLE  4-1.  A1,  A2,  u1  for  e 

e  e  e 


0.1,  0.04,  0.01,  0.001  | 

I 

i 


e 

A1 

5? 

e 

£ 

e 

0.1 

0.77230 

4.5742 

2.3462 

0.04 

0.36160 

1.3715 

2.4782 

0.01 

0.09650 

0.38677 

2.2153 

0.001 

0.01069 

0.03997 

2.2538 

a 

It  is  clear  that  most  of  the  features  described  in  Section  3  are 

exhibited  in  these  plots.  First,  observe  the  attenuation  of  k  *  1,2  as 

e 

e  +  0  on  Qq.  Second,  note  the  oscillatory  behavior  of  ^  as  e  0  on  Qj. 

e 

Furthermore,  one  may  add  that  the  first  eigenvalue-eigenvector  pairs  are 

computed  more  accurately  than  their  last  counterparts  [13].  Consequently,  if 

the  eigenvalues  of  A*1  are  ordered  ascendingly,  is  pushed  higher  and 

e  e  e 

higher  as  e  +  0  and  hence  computed  less  and  less  accurately.  Moreover,  since 

h  is  fixed,  the  oscillatory  behavior  of  \|A  would  not  be  captured  by  this 

e 

approximation,  unless  h  is  made  smaller  and  hence  increasing  the  order  of  the 
matrix  A*1. 


:1 
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SECTION  5 


CONCLUSION 

The  eigenvalue  problem  of  a  class  of  stiff  operators  has  been  analyzed  in 
this  paper,  via  a  general  formulation  using  bilinear  forms,  to  avoid  the  com¬ 
plexity  of  explictly  keeping  track  of  the  various  boundary  and  interface 
conditions.  First,  the  intuitive  idea  that  the  eigenvalues  of  stiff  operators 
are  of  different  order  of  magnitude  as  functions  of  the  parameter  e,  is 
rigorously  verified.  Second,  a  terminology  is  introduced  to  describe  the  con¬ 
vergence  of  the  eigenvector  as  e  ♦  0,  i.e.,  flattening,  alternation  and  oscil¬ 
lation.  This  analysis  is  of  paramount  importance,  because  it  will  yield  in¬ 
sight  into  how  to  approximate  boundary  value  problems  involving  stiff 
operators. 
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ABSTRACT 


In  this  paper,  asymptotic  approximations  for  a  class  of  linear  stiff 
systems  are  constructed,  using  the  weak  limits  of  the  eigenvectors  of  the 
associated  stiff  operators. 
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SECTION  1 


§ 

vw 


INTRODUCTION 

A  standard  approach  in  solving  many  boundary  value  problems  of  mathemat¬ 
ical  physics  is  to  employ  eigenvector  expansions.  When  such  problems  contain 
one  (or  several)  small  parameter  e,  asymptotic  expansions  in  (fractional) 
powers  of  e  can  often  be  a  powerful  approximating  process. 

In  a  previous  paper  [5 ] * ,  the  spectral  decomposition  of  a  class  of  stiff 
operators  was  investigated.  The  asymptotic  properties  of  their  eigenvalues 
and  the  corresponding  eigenvectors  were  analyzed.  More  importantly,  it  was 
discovered  that  some  of  the  aforementioned  eigenvectors  are  not  analytic 
functions  of  e.  This  fact  is  intimately  related  to  their  convergence  in  a 
Hilbert  space  with  a  weaker  topology.  Despite  this  undesirable  property,  it 
is  possible  to  derive  asymptotic  approximations  of  the  solutions  of  boundary 
value  problems  involving  stiff  operators. 

In  the  sequel,  the  results  of  [5]  are  employed  to  investigate  the  beha¬ 
vior  of  the  solutions  of  the  following  three  formal  equations 

t 

Ae  ye  -  f  (1-1) 

3  ye 

“ —  +  Ae  ye  -  f  ,  ye(0)  -  h  (1-2) 

3t 


♦References  are  indicated  by  numbers  in  square  brackets,  the  list  appears  at 
the  end  of  this  Appendix. 
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(1-3) 


32ye  3Ye 

-  +  Ae  ye  =  f  .  ye(0)  *  h  ,  -  (0)  =  g 

3t2  3t 

where  ^  is  a  stiff  operator  of  the  class  studied  in  [5). 

The  occurrence  of  ( 1—1 )— ( 1—3 )  is  very  frequent  in  mathematical  models 
of  distributed  physical  processes  such  as  nuclear  reactors,  heat  exchangers, 
chemical  reactors,  fluid  systems,  vibration  systems,  steel  and  glass  proces¬ 
ses,  etc.  Thus,  it  is  important  to  focus  on  them.  The  major  thrust  of  this 
paper  is  to  derive  asymptotic  approximations  of  the  solutions  ye  of  (1-1)- 
1-3).  For  elliptic  problems  (i.e.,  (1-1)),  by  appropriately  modifying  the 
weak  limits  of  the  eigenvectors ,  one  may  be  able  to  compute  an  asymptotic 
expansion  of  any  order  for  ye.  In  so  doing,  the  formal  results  obtained  in 
[3]  are  complemented.  For  revolution  problems  (i.e.,  (l-2)-(l-3)) ,  the  con¬ 
cept  of  weak  solutions  [2], [4]  is  used  to  define  zeroth  order  approximations 
of  the  solutions  of  ( 1—2 )— ( 1—3 ) ,  which  depend  only  upon  the  weak  limits  of 
the  eigenvectors. 

This  paper  is  organized  as  follows.  In  Section  2,  the  main  results  of 
(5]  are  briefly  reviewed.  In  Section  3,  the  solution  of  (l-l)  is  derived. 

In  Section  4,  the  convergence  of  the  solution  of  (1-2)  is  investigated  and  an 
asymptotic  approximation  is  constructed  for  it.  In  Section  5,  an  analysis 
similar  to  that  of  Section  4  for  a  class  of  hyperbolic  problems  (i.e.,  (1-3)) 
is  undertaken.  In  the  last  section,  some  concluding  remarks  are  presented. 


<> 
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SECTION  2 


REVIEW  OF  PREVIOUS  RESULTS 


Let  V,  H  be  two  given  real  Hilbert  spaces  such  that  V  is  dense  in  H  and 


Al)  the  injection  of  V  into  H  is  compact. 


Let  V*  denote  the  dual  space  of  V.  After  identification  of  H  with  H*,  one  has 


A2)  VCHCV*. 


Let  i  *  0,1  be  two  forms  on  V  such  that  the  following  assumptions 


hold: 


A3)  ai(<|>l'l')  is  bilinear,  symmetric  on  V 


A4)  a ( 4» , ^ )  is  contiuous  on  V,  i.e.,  there  exists  8i  such  that 


ai(<|),<|>)  <  8i  8$ By  “’Hv,  V<J»eV,  Vi|»eV 


A5)  >  Oipi(^)2,  where  a±  >  0  and  pi(*)  is  continuous 

semi-norm  on  V. 


A6)  po(<J>)  +  Pl($)  is  a  norm  equivalent  to  l$iy 


A7)  a^ ( , <}> )  -  0  on  VjCV,  where  V^  is  an  infinite-dimensional 
subspace  of  V,  i«0,l. 


A8)  If  \|/  H-Lo(^)  is  a  continuous  linear  form  on  V,  null  on  Vo, 
there  exists  $eV  (module  Vq)  such  that 


«()<♦.*>  ’  I*0( *> •  **eV 


Let  ae(ij),\|))  be  defined  as 


ae( ♦ » ♦)  *  a0( +  e  ai(<M) 


(2-1) 


Remark  2. I : 


Clearly,  the  form  ac(4,<|;)  is  bilinear,  symmetric,  coercive  and  bounded 
on  V.  Hence,  it  defines  uniquely  a  selfadjoint  operator  Ag 

Ae  e  L(V;V*)  , 

aeU,<|>)  “  <Ae  <M>  .  for  all  4» , £  V 

where  <*,*>  denotes  the  duality  pairing  between  V  and  its  dual  V*. 


From  the  above  remark,  one  concludes  that  the  spectrum  of  Ae  is  a  subset 
of  ]R+,  consisting  only  of  the  point  spectrum  [1]. 

00 

The  eigenvalue  problem  for  A£  is,  then,  to  seek  {Yk»Xk)  G  x  V 

g  e  k=l 

such  that 


A  Xk  =  Yk  Xk  (2-2) 

e  e  e  e 

or  equivalently 

a  (xk»$)  “  Yk(xk»,fr)  »  for  all  y  e  V  (2-3) 

e  e  e  e 


Now  the  properties  of  (Yk»Xk)  are  summarized  in: 

e  e  k**l 


Theorem  2» 1 : 

00  00 
Let  (yk)  be  the  eigenvalues  of  Ae  and  (xk}  the  corresponding  nor- 
e  k»l  e  k»l 

malized  system  of  eigenvectors.  Then,  given  a  sequence  of  e  converging  to 
00  00  00 

zero,  (Yk,Xk)  can  be  decomposed  into  {Xk,$k}  and  (yk,^k)  such  that, 
e  e  k»l  e  g  k*l  g  g  k=l 

for  each  k, 
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1)  Ak  ■*  0  linearly  in  £,  ^  +  ^  weakly  in  V 
e  e 


2)  pk  yk  >  o,  ijik  ■*■  ij)k  weakly  in  H 

0  e 

OO  CO 

where  {<{>k}  and  {tfik}  satisfy 
k=l  k=l 


a  (<J>k,x)  =  Ak(<J>k,X),  *k  e  V  C  V,  for  all  x  e  V  (2-4) 

110  0 


a  (4»k»x)  =  yk(\j»k,x) ,  i(»k  e  H  C  H,  for  all  x  e  v  (2-5) 

0  0  1  1 


Moreover,  the  following  estimates  hold  for  e  >  0: 

3)  l*k«  <  C  (2-6) 

e  y  1 

4)  FT  l«l>k«  <  C  (2-7) 

e  V  2 


where  C[,  C2  are  two  positive  constants,  independent  of  e. 


Proof:  See  [5]. 


Remark  2.2: 


The  subspace  Hj  of  H  in  (2-5)  is  uniquely  determined  as  the  span  of 

"  . 

k*l 


Remark  2.3: 

OO  A  OO 

The  weak  limits  {d>k}  and  {ij>k}  form  an  orthonormal  system  in  H. 

k-1  k-1 


SECTION  3 


ELLIPTIC  BOUNDARY  VALUE  PROBLEMS 


In  this  section,  an  asymptotic  expansion  of  the  solution  ye  of 


Ae  ye  =  f  ,  f  e  H  (3-1) 

or  equivalently 

a0 C y e » 4> )  +  e  ai(ye»4>)  “  (£,♦)  .  e  V  ,  for  all  $  e  V  (3-2) 

is  desired.  This  problem  is  studied  in  a  more  general  context  (namely  without 
assumption  Al)  in  [3]. 

The  modifications  of  the  weak  limits  of  the  eigenvectors  are  now  consid¬ 
ered.  Since  ^  e  H  ,  one  may  add  to  it  a  function 
e  1 

-  &  +  e  +  c2  +  ...  (3-3) 

e  0  1  2 


where 

c  H  such  that  ^  e  V  (3—4 ) 

0  0 

(JeV  ,  l  -  1,2 .  (3-5) 

This  amounts  to  solving  the  following  boundary  value  problem,  using  the  iter¬ 
ative  process  of  [3]: 

a  +  e  a  ('Ji^+C'Sx)  ■  P^(^»x)  »  for  all  x  e  v  (3-6) 

0  e  1  e  0 
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Consequently,  ,{£k}  satisfy 
*  £=0 


a  (Ck,x)  =  0  ,  ^  +  C'1  e  V  ,  for  all  x  e  V 

10  0  0 


(3-7) 


3  (€k,x)  =  “  a  (£k  ,x)  ,  for  all  x  e  V 

0  £  1  £-1 


3  (ek»x)  =  0  ,  for  all  X  e  V 

1  £  0 


(3-8) 


£  =  1,2,... 


k  =  1,2,... 


Remark  3.1: 


It  is  worthy  of  mention  that  the  iterative  process  described  by  (3-7)- 

(3-8)  appears  to  average  the  oscillatory  behavior  of  i|/k  [5], 

e 

Similarly,  one  adds  to  $k 


0k  »  e  ek  ■+  e2  Qk  +  . . . 
e  1  2 


(3-9) 


where 

0k  e  V  ,  £  *  1,2,... 

*  , 

The  function  0k  is  chosen  such  that 
e 

a  (4>k  +  0k,x)  +  e  a  ( 4>k, x)  *  e  lk(<i>k»x)  »  for  all  x  c  v 
0  e  1  1 


from  which  one  concludes  (using  the  fact  that  <j>k  e  Vo)  that  {0kJ 

£  £*l,k=l 

satisfy 


a  (0k,X>  =  a  ( 4»k, x)  .  for  all  x  e  V 
0  1  1 


a  (9k,x)  3  0  *  for  all  X  c  V 

1  1  0 


a  Ok,x)  =  -  a  <6^  ,x)  ,  for  all  x  e  v] 

04  14-1 

ai(8k,x)  =  0  »  for  all  x  e  Vq 

4 


9  =  2,3,... 
k  *  1,2,... 


Remark  3.2: 


The  zeroth  term  in  (3-9)  is  zero  because  $k  e  V. 


(3-10) 


Remark  3.3: 

The  iterative  process  presented  above  appears  to  average  the  flattening 

00 

(and  attenuation)  of  {$k}  [5]. 

e  k=l 

Before  proceeding  further,  consider  the  following  example  to  examine  what 
(3-7 )-(3-8) ,  and  (3-10)  yield: 

Example  3.1: 

Let  H  *  L 2(0),  V  *  Hl(Sl)  where  SI  ■  fin  U  ^1  ^  S,  T  ■  Tq  U  I*i  as  indicated 

0 

in  Fig.  3-1.  Let 

n  34>  3<J» 

ai( ♦,♦)  ■  E  J  — -  - —  dx 

j*l  flt  3xj  3xj 


B-10 


in'  mmwi 


Figure  3-1.  Examples  of  Interfaced  Sets 


then  (3—7 )— (3— 8) ,  and 

(3-10) 

become 

£k  -  o 

on  0 

00 

0 

-Hk  *  0 

on  0 

01 

1 

Sk 

-  0 

,  5k 

* 

01 

r. 

01 

s  0 

1 
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EJP16X5W  kVA.%.  klv 


/»' 


-A  0k-  =  0  on  ft 

£0  0 


£0 


36 


=  0  , 


£0 


3v 


36 


£-1,1 


3v 


-A  8k  -  0  on  ft 
£1  1 


8k 

=  0  , 

8k 

_ 

ek 

£1 

F1 

£1 

S 

£0 

£  *  2,3,... 
k  =  1,2,... 


Remark  3.4: 

The  computation  of  these  modifications  is  greatly  simplified  since  the 
dependence  on  the  parameter  e  is  eliminated. 

Now  attention  is  focused  into  obtaining  an  asymptotic  expansion  of  the 
solution  of  the  boundary  value  problem  (3-1)  using  the  modified  weak  limits 
of  the  eigenvectors  of  Ag. 

Theorem  3.1: 

l 

For  sufficiently  small  e ,  the  solution  of  (3-1)  is  given  by 

*  ck  •  dk 

ye  -  S  —  (*k+5k)  +  E  -  (*k+0k)  (3-11) 

k-1  yk  e  k-1  eXk  e 
0  1 
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where 


{pk^k}  and  {A*,**}  satisfy  (2—4 )— ( 2— 5 ) 
0  k“l  l  k=l 


{^k}  and  {a*1}  are  given  by  (3-3)  and  (3-9) 
e  k=l  e  k*l 


{cM  and  (dk)  are  the  Fourier  coefficients  of  f,  i.e., 
k-1  k=l 


ck  «  (f^k)H 


(3-12) 


dk  *  (f,^)H 


(3-13) 


Proof: 


The  eigenvectors  (ik)  and  {i|>k}  form  a  complete  orthonormal^  basis 

k-1  k-1  - 


of  H  and  <|>k  +  0k  e  V,  i)ik  +  e  V  by  construction.  Now  it  is  straightforward 
e  e 


to  verify  that  (3-11)  is  the  unique  asymptotic  expansion  of  ye. 


Remark  3.5: 


Truncate  £k  and  0k  to  the  pth  term  and  denote  these  truncated  series, 


respectively,  by  (k,p  and  0k»P.  Definte  eP  by 
e  e  c 


eP  ■  ye  -  yp 

e  e 


where 


*  ck  •  dk 

yP  -  I  —  +  (k,p)  +  E  -  (4>k  +  0k, P)  . 

e  w  uk  e  k-l  exk  e 

0  l 


^Note  that  must  be  renormalized  in  H, 

k-1 


t  J%E  Jfli  i  TV  %J%1 


Then  it  can  be  shown  ([3],  p.  13)  that 


JePI  <  C  eP 
e 

where  C  is  a  constant  independent  of  e. 


Example  3.2  (Example  3.1  Cont.) 

In  this  case,  (3— 11)— (3— 13)  become 


“  ck  *  dk 

y  =  £  —  (\{/k+£k  )  +  £  -  gk 

k^l  pk  0  eO  k=l  elk  fO 

0  1 


•  ck  “  dk 

f  *  I  -  +  E  -  (^k+gk  ) 

el  uk  el  it»i  exk  1  el 

0  1 


where 


ck  *  (f  ,\J;k) 


L2(«) 


<f 

0  °L2(n0) 


dk  «  (f,i|»k)  -  (f  ,^k) 

L2(n)  1  1  L^Cflj) 


i,k  ,  0k  are  computed  in  Example  3.1. 
e  e 


B-15 


SECTION  4 


PARABOLIC  BOUNDARY  VALUE  PROBLEMS 


In  this  section,  an  evolution  problem  of  parabolic  type  is  considered. 
Hence,  let  the  variable  t  denote  time.  It  is  assumed  that  t  e  (0,T),  T  <  »  , 

and  that  all  the  assumptions  made  in  Section  2  hold.  Let  L2(0,T;V),  L2(0,T;H) 
L2(0,T;V*)  denote  the  Hilbert  spaces  of  Lebesgue  square  integrable  functions 
with  values  in  V,  H,  V*',  respectively.  Let  prime  denote  the  distributional 
derivative  with  respect  to  time  [4].  In  the  sequel,  the  following  parabolic 
boundary  value  problem  Is  analyzed 


(y',<j>)  +  a  (y  ,$)  +  e  a  (y  ,$)  -  (f,<j>)  ,  for  all  <(>  e  V  (4-1) 

e  0  e  l  e 


y  (0)  =»  h  ,  h  given  in  H 
e 


(4-2) 


y  e  L2(0,T;V)  ,  f  e  L2(0,T;H)  . 

e 


(4-3) 


Under  the  assumptions  made,  problem  (A— 1 )— (4— 3 )  admits  a  unique  solution 
ye  e  L2(0,T;V)  (4j.  Using  the  results  derived  in  Section  2,  the  convergence 
of  ye  as  e  0  is  studied.  Then  an  asymptotic  approximation  of  ye  is  con¬ 
structed  and  an  asymptotic  error  estimate  is  derived. 


4. 1  CONVERGENCE  OF  yr  as  e  ♦  0 

at  at 

As  in  Section  2,  let  Uk,*k}  and  {pk,\J»k)  be  the  exact  eigenvalue- 

c  e  k-1  e  e  k-1 

eigenvector  pairs  of  Ae  with  the  eigenvectors  normalized  to  one  in  H. 


(4-4) 


1c  ■  1,2,... 

whose  solutions  are 


-ykf 

1 

t 

-jJk(t—r) 

ck(t)  *  e  e 

(h,ij^) 

+  / 

e  e  ( f ,  i|ik)  dr 

e 

e 

0 

e 

-Xkt 

t 

-Xk(t-T) 

dk(t)  -  e  e 

(h,*k) 

+  / 

e  e  (f,$k)  dT 

e  eg  c 


Using  Theorem  2.1,  one  concludes  with  no  difficulty  that 

ye  y  weakly  in  1.2(0, T;H) 


where 


y  *  I  ck  <|>k  +  S  dk  $k 
k*l  k=l 


{*k}  {^,k}  are  (-he  Units  of  (4>k)  (respectively,  {^k}  ) 

k=l  k=l  e  k=l  e  k=l 

in  V,  (respectively,  in  H)  given  by  (2-4)  (respectively,  2-5) 


-Ukt  t  —  pk(  t— t  ) 

ck(t)  »  e  0  (h,4»k)  +  /  e  0  (f,i|»k)  dr 

0 


dk(t)  -  (h,$k)  +  /  (f,$k)  dT  . 

0 


The  preceding  discussion  is  summarized  in: 

Theorem  4,1: 

Let  ye  denote  the  solution  of  (4-l)-(4-3).  Given  a  sequence  of  e  con¬ 
verging  to  zero, 

ye  ♦  y  weakly  in  L2(0,T;H)  (4-8) 

where  y  is  given  by  (4-7).  Moreover, 

>/r  ly  I  <  C  (4-9) 

.  E  L2(0,T;V) 

Where  C  is  a  constant  independent  of  e. 


Proof:  ’ 

Use  the  eigenvalue-eigenvector  pairs  of  the  operator  Ae  and  their  proper- 

\ 

ties,  as  described  by  theorem  2.1,  to  obtain  (4-8).'  The  estimate  (4-9)  is 
then  readily  derived  by  employing  (2-6)-(2-7). 
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4.2  ASYMPTOTIC  APPROXIMATION  OF  y 


It  should  be  noted  that  the  method  by  which  stiff  elliptic  boundary  value 
problems  were  solved  in  Section  3  does  not  yield  an  iterative  process  for  evo¬ 
lution  problems  in  general.  Therefore,  one  would  be  content  to  obtain  a  “weak" 
approximation  of  the  solution  of  (4—1 >—(4—3 )  using  only  the  weak  limits  of  the 
eigenvectors  of  Ag. 

In  the  sequel,  an  approximation  of  the  solution  of  the  following  boundary 
value  problem  is  derived: 

3yeO 

— —  +  Ao  ye0  ■  fo  on  Qo 


3yel 

— ■  +  e  Ai  yel  *  fx  on  Qi 


yeQ  *  0  on  Iq  ,  yei  *  0  on  Z\ 


(4-11) 


yeo  *  yel 


3ye0  3yel 

3%  3vAl 


on  R 


(4-12) 


ye(0)  *  h  on  8 


(4-13) 


where 


Qi  -  «i  x  (0,T)  ,  i  -  0,1 
li  -  rt  x  (0,T)  ,  i  -  0,1 
R  “  S  x  (0,T)  . 
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a^  satisfy 

ij 


1)  a k  e  Ckflk) 
ij 


2)  ak  =  ak 

ij  ji 

n  n 


3)  £  E  aM  5  >  <»lc  i  5  ,  «  >  0 

i-1  j,i  ij  i  j  i-1  i  k 

for  all  5  £  Rn  ,  5*0  . 


The  chief  reason  for  this  digression  is  that  one  can  explicitly  specify 
the  regularity  conditions  of  the  functions  involved  in  the  construction  of  the 
approximation  of  ye.  However,  the  concepts  involve  herein  are  equally  appli¬ 
cable  to  the  general  case. 

Now,  let  the  zeroth  order  approximation  be  denoted  by  y°  ■  (y®,y°  )  and 

e  0  el 

defined  as  follows: 


\ 


3y0 

-  +  A  y0  »  f  on  Q 

at  ooo  o 


y0  *  0  on  E 
0  0 


(4-14) 


*yr 


3v 


0  on  R 


AO 


y0(0)  -  h  on  0 

0  0  0 


'<1 


w 


-  +  c  A  yO  =  f  on  Q 

3t  1  el  1  l 


yy  =  0  on  £ 
el  1 


(4-15) 


yu  on  R 
0 


y°  (0) 

el 


h  on  ft 
1  0 


Remark  4. 1 : 


The  solution  of  y^ 'of  (4-14)  is  regular.  Actually,  y^  e  L2(0,T;H*(ft^; T  ) ) 


CL2(Q  ).  Hence  y^  e  L2fO,T;H^/2(S) )  CL2(R).  Problem  (4-15)  is  a  nonhorao- 
0  0  R 


geneous  boundary  value  problem.  Consequently,  y^  has  meaning  in  a  weak  sense 

el 


using  transposition  [4], 


Since  the  zeroth  order  approximation  y0  is  weak,  one  rewrites  (4-14)- 

e 


(4-15)  in  the  proper  form,  using  transposition  as  follows.  Let 


>1  ■  <  Xi  :  Xi  e  h2.1(Q!)  , 


o 

X 

•  0  >  : — 

-  0  ,  X0(T)  -  0 

ZQ  3vAo 

R 

-  o  ,  X! 

-  0  ,  XL(T)  -  0 

z, 

R 

1 

Now  consider  the  following  isomorphisms 


X 


0 


*0 


3x0 

-  +  Ao  xO  from  *0  to  L2(Qo) 


*1 


£1 

I— 


3*1 

-  —  *  c  A, 


*1  from  to  L2(Qi) 


By  transposition,  one  concludes  that  Xg  I—*  M®(Xg),  being  a  continuous  linear 
form  on  *g  (endowed  with  the  topology  induced  by  H2»l(Qg)),  there  exists  a 
unique  y®  e  L2(Qq)  such  that 


I 

Qo 


dQo  ■  M°(Xq) 


for  all  Xg  e  *g 


(4-16) 


and 


m0  | — .  y0  is  a  continuous  linear  mapping  of 
0 


0*  L2(Qg)  (cf.  Remark  4.2). 
TO 


Similarly,  x  I — *  M^(x  )  being  a  continuous  linear  form  on  4>  (endowed  with 
lei  1 


the  topology  of  H2**(Q  )),  there  exists  a  unique  y®  e  L2(Q  )  such  that 


el 


Xl 

f  y0  (  -  -  +  eA  x  1  dQ  ■  M^(x  )  »  for  all  X  G  *  (4-17) 

ql  el  \  3t  11/0  e  l  11 


and 


M*  1 — ►  y0  is  a  continuous  linear  mapping  of  #*  I — •  L2(Q  ) 
e  el  11 
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Select  )  and  mUx  )  as 

0  e  1 


M°(X  > 
0 

MJ(X  ) 
e  1 


/  fnX  dQ  +  /  h  x  (x.O)dft 
0  0  0  o  0  0  C 

*0  0 


/  f  x  <*Q  +  /  1>  X  (x,0)dft  - 
Ql  1  1  0  1  1  1 


e  /  yO  -  dR 

R  0 


Remark  4.2: 

Since  Cf*(Q  )  (L  ♦  ,  the  dual  of  ♦  is  not  a  space  of  distributions. 

0  0  0  0 

Therefore,  the  introduction  of  0  is  required  to  interpret  duality  [4], 

ro 

It  can  be  easily  verified  that  the  solutions  of  (4—1 4 )— (4— 1 5 )  and  (4-16)- 
(4-17)  are  identical  [2]. 

Now  an  error  estimate  between  the  exact  solution  ye  and  its  approximation 

y0  is  derived, 
e 


Theorem  4.2: 

Let  ye  be  the  solution  of  (4— 15)— (4— 16)  and  y®  the  solution  of  (4-16)- 

e 

(4-17).  Then,  for  sufficiently  small  e,  one  has 

<ye  ~  y°»  <  C  e1/2  (4-18) 

>£  L2(Q) 

where 

Q  =»  Qq  x  Ql 

and  C  is  a  positive  constant.  Independent  of  c. 

Proof: 

First,  rewrite  (4— 10)— (4— 13)  in  the  weak  form,  i.e., 
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(4-19) 


/  9*o  \ 

(  3t  +  Vo) 


dQ  *  /  f  X  <*Q  +  J  h  X  (x,0)dfi 
0  0  0  0  0  0  0 

W0  0 


/  j  tv 

-  Y  dR  ,  for  all  x  e  * 

R  3%  »  0  0 


/  3X1  \ 

( - +  e  A  X  dQ  s  |  f  X  <i Q  +  |  l>  X  (x,0)dfl 

\  3t  1  ll  1  Q,  1  1  1  0  1  1  1 


(4-20) 


-  e  /  y  -  dR  ,  for  all  x  e  $ 

R  £1  3vA!  •  l  l 


Subtract  (4-16)  from  (4—19)  and  (4-17)  from  (4-20)  to  obtain 


/  n\  /  3X°  \  , 

| y  *2  i "  r” +  dQrt " e  / ; —  dR , 

leO  0  /  1  3t  00/  0  3vAl  0 


for  all  x  e  * 

0  0 


(4-21) 


J  (y  ~  y®\( "  +  e  A  x  jdQ  ■  e  /  (y  -  y®\  -  dR  ,  for  all  x  e  * 

i  V*1  “A  u  “/  1  R  V*0  °/  »*»i  1  1 


(4-22) 


where  the  interface  condition  (4-12)  is  used. 


Now  consider  the  following  equations: 


3Xe° 

-  +  Ax  -  y  -  r  on  Q 

3t  0  cO  eO  0  0 


X  *  0  on  E 
eO  0 


(4-23) 


0  on  R 


X  (x,T)  -  0 
eO 


3Xel 

-  +  e 

A  X 

3t 

0  el 

X 

11 

o 

on  I 

el 

l 

X 

II 

o 

on  R 

el 

X  (x,T) 

=  0 

y  -  y®  on  Q 


(4-24) 


el 


Since  the  coefficients  of  A^,  i  *  0,1  are  assumed  to  be  sufficiently  regular 


and  y  -  y°  e  L2(0,T;L2(fl  )),  y  -  y0  e  L2(0,T;L2(n  ))  ,  one  concludes  that 
e0  0  0  el  el  1 


X  e  *  ,  i  -  0,1  . 

ei  i 


Let  X  *  X  in  (4-21)  and  x  *  X  in  (4-22),  and  use  (4-23)-(4-24)  to  obtain 
0  eO  1  el 


•y  _ 

■  yJ" 

<  Cj  e1/2 

eO 

0  l2<V 

»y  ,  * 

-  y°  « 

<  C2  e1/2 

el 

L2(Qj) 

C^,  are  some  positive 

constants , 

independent  of  e.  Hence, 

xT 

i 

o 

i 

> 

1 

-  y0< 

+  ly  1  <  C  e1/2 

e  e  l2(q)  eu 

0 

el  l2(q1> 

Remark  4,3: 


It  should  be  noted  that  (4-18)  holds  for  e  small,  but  strictly  posi¬ 
tive,  because  of  the  heavy  reliance  upon  the  fact  that  the  solution  of  (4-24) 


B-25 


3 


belongs  to  It  is  easy  to  see  that  if  one  formalLy  sets  e  »  0  in  the  first 

equation  of  (4-24),  its  right-hand  side  would  be  in  L,2(o,T;  1,2(0  ))  but  x 

1  el 

for  e  =  0  would  not  be  in  4  . 

I 

Recall  that  it  was  shown  in  Section  2  that  the  eigenvalue-eigenvector 

00  00 

pairs  of  A  ,  i.e.,  {y^,xk)  are  decomposable  into  two  groups  and 

e  e  e  k=l  e  e  k=l 

{pk^k} 

e  e  k“l 


The  normalized  weak  limits  (in  hI(Q))  of  satisfy 


e  k-1 


»  0  on  ft 
0  0 

A  $k  «  \k  *k  on  q 
1111  1 


<j>k 

1 


0,  ^ 


(4-25) 


k  *  1 ,2, . . • 


The  normalized  weak  limits  (in  1,2(0))  of  {\|)^}  obey 

e  k=«l 


A  \(»k  «  yk  ^k  on  0 
0  0  0  0  0 

j,k  ■  0  'on  0 
0  1 


0 


A 


o. 


3v 


AO 


(4-26) 


k  -  1,2,... 


Then  the  solution  of  (4-14)-(4-15)  may  be  represented  by 
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V 


yO  =  E  ck(t)r|ik 
0  (c=l  0 


■\ 


00  f 

ck  e  L2(0,T),  E  /  |ck(t)|2  dt  < 
k=l  0 


yO  =  E  dk(t)$k 
el  k=i  e  1 


CO  T 

dk  e  l2(0,T),  E  /  |dk(t)|2  dt  <  « 
e  k-1  0  e 


In  order  to  determine  ck(t),  dk(t),  let 

e 


X  (x,t)  =  0(t)<pk(x),  0  e  Cl([0,T]),  0(T)  -  0  (so  that  x  e  #  ) 
0  0  0  0 


X^(x,t)  -  v(t)$k(x),  v  e  C^(  [0,T] ) ,  v(T)  -  0  (so  that  e  M 


in  (4— 16)— (4— 17)  to  get 


T 

/  ck 
0 


d0  .  \  fT 
—  +  yk  I  dt  ■  /  (f  ,i|>k)  0  dt 
dt  0/  o  0  0 


(h  ,<|>k)0(O) 
0  0 


T  /  dv  \  T 

/  dk - +  Xk  e  v  dt  -  /  (f  ,$k)v  dt  -  (h  ,$k)v(0) 

0  e  \  dt  l  /  0  1  1 


I”! 


T  /  •  "yj 

-  e  /  I  E  cA  /  <|»*  — —  dS  J  v  dt  , 
0  U-l  S  3vA[ 


which  are  equivalent  to 
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dck  ' 

—  +  ukck  =•  (f  ,<|»k) 
dc  0  0  1 


ck(0)  =  (h  ,ij/k) 
0  0 


d  dr 


-  +  Ak  e  dk 

dt  1  e 


<f  ,^) 


3*, 


e  E  c ^  f  ^  -  dS 

4-1  S  0  3vAj 


dk(0)  =  (h  ,*k) 
e  1  1 


Remark  4.4: 


Using  the  terminology  of  singular  perturbation,  a  two  time-scale  decom¬ 


position  is  achieved.  Namely,  y0  is  “fast”  and  yQ  is  predominantly  "slow", 

0  el 


B-28 


SECTION  5 


m 


HYPERBOLIC  80UNDARY  VALUE  PROBLEMS 


In  this  section,  an  evolution  of  hyperbolic  type  is  investigated,  namely, 


the  following  boundary  value  problem: 


(y",x)  +  a  (y  .x)  +  £  a  (y  »x)  ■  (f»x)»  for  all  x  e  V 
e  0  e  1  e 


(5-1) 


y  (0)  -  h,  h  given  in  V 
e 


(5-2) 


y'(0)  -  g,  g  given  in  H 
e 


(5-3) 


y  e  L2(0,T;V),  y*  c  L2(0,T;H),  f  e  L2(0,T;H)  . 

e  e 


(5-4) 


The  present  analysis  will  be  parallel  to  that  of  Section  4,  Hence,  the 


convergence  of  y  as  e  ♦  0  is  studied.  Then  a  zeroth  order  approximation  of 
e 


y  is  constructed,  using  the  weak  limits  of  the  eigenvectors  of  A  (i.e., 
e  e 


the  operator  associated  with  the  bilinear  form  a  ($,<)/)  ■  a  (4>,^)  +  c  a  (<J>,iJ>) 

e  0  1 


(cf.  Section  2).  As  previously  indicated,  (5-l)-(5-4)  have  to  be  specialized 


to  second  order  operators  A  ,  so  that  one  can  specify  precisely  what  is  re- 

e 


quired  for  the  present  analysis  to  hold.  Therefore,  the  following  problem  is 


considered  in  the  sequel: 


-  —  +  A  y  ■  f  on  Q 
3t2  0  e0  0  0 


-  +  e  A  y  *  f  on  Q 

3t2  1  el  l  1 


(5-5) 


(5-6) 


y  =  0  on  E  ,  y  =■  0  on  E 
eO  0  el  l 


y  "  y  , 

eO  el 


3yeO  3yel 


3 


(5-7) 


y  (0)  -  h  ,  h  -  (h  ,h  )  ,  h  e  HKO  )  ,  i  =  0,1 

e  0  1  i  0  i 


W  Q  / 


-  (0)  -  g  ,  g  e  L2(0) 

3t 


(5-9) 


y  e  L?(0,T;  Hl(fl))  ,  -  e  L2(0,T;  L2(Q)) 

e  0  3t 


(5-10) 


f  e  L2(0,T;  L2(0) ) 


(5-11) 


where  A^,  1-0,1  satisfy  the  conditions  of  subsection  4.2 


5. 1  CONVERGENCE  OF  y  AS  e  ♦  0 

e 


Theorem  5.1: 


Let  y  be  the  solution  of  (5-5)- (5-11 ). 
e 


Then,  given  a  sequence  of  e 


converlng  to  aero, 


y  ♦  weakly  in  L2(0,T;L2(G>) 
e 


(5-12) 


3yC  3y 

—  ♦  —  weakly  in  L2(0,T;L2(fl) ) 
3t  3t 


(5-13) 


m 


Moreover , 


y/z  II y  1  <  C  (5-14) 

£  L2(0,T;  H*(R)) 

0 


C  is  a  constant  Independent  of  e. 


Proof: 


Let  {Ak,$k}  and  be  the  exact  eigenvalue-eigenvector  pairs 

e  c  k-J.  e  e  k-1 

of  A  ,  with  the  eigenvectors  normalized  (in  L?(R))  for  fixed  e.  Using  a 
e 

finite  dimensional  approximation  of  y  such  as 

e 


m  m 

y®  =  £  ck  +  £  dk 

e  k-1  e  e  k-1  e  e 


is  is  shown  in  [2]  that 


ym  +  y  strongly  in  L2(0,T;  HUR))  as  m  ® 
e  e  0 

3y"  3ye 

-  +  -  strongly  in  L2(0,T;L2(R))  as  m  +  +  ■ 

3t  3t 

I 

where  , 


y  -  Z  c^  ft  +  £  d^  ft  (5-15) 

e  k»l  e  e  k-1  e  e 


(cty  ,  (dM  satisfy  the  following  set  of  ordinary  differential  equations 
c  k-1  e  k-1 


d2ck 

- —  +  pk  ck  a  (f,<|ik) 

dt2  e  e  e 


ck(0)  -  (h,*k> 
e  e 


—  (0)  -  (g,<|)k) 
dt  e 


d2dk 

- -  +  Xk  dk  «  (f,*k) 

dt2  e  e  e 


dk<0)  -  (h,<>k) 
e  e 


—  (0)  *  (g,$k) 
dt  e 


whose  solutions  are  given  respectively  by 


ck(t) 

e 


__  ^.♦g)  t  sin ‘y^Ct-T) 

(h,ij>k)  cos  +  - -  sin  t  +  /  - - - -  (f,\|»k)dx 

e  ■  * e  ▼  e  0  E 


(5-16) 


dk(t) 

e 


(g.^> 


i,$k)  cos  -/I^t  +  — — sin 
e  '  c 


_  t  8in  V*«> 

/xkt  +  / - (f,$k)dx 

*  c  o  e 


(5-17) 


Now,  let  e  be  a  sequence  converging  to  zero.  The  convergence  property  (5-12) 
and  the  estimate  (5-14)  are  easily  deduced  from  the  results  of  Theorem  2.1, 
where  y  is  written  as 


00  00 

y  =  l  ck  +  i  <jk  $k  (5-18) 

k=l  k=l 


where  (ck)  ,  {dk} 

k-1  k=l 


satisfy  (5— 16)— (5— 17)  after  letting  e  +  0, 


i«6«  ) 


T 

dk(t)  -  (h,$k)  +  (g,<J>k)t  +  /  (t-t)(f ,4>k)dT  . 

0 


Differentiate  (5-15)  with  respect  to  time  and  take  limit  as  e  +  0  using  the 
results  of  Theorem  2.1  to  get  (5-13). 

5.2  ASYMPTOTIC  APPROXIMATION  OF  y 

e 

In  the  sequel,  the  zeroth  order  approximation  y®  of  y  is  constructed 

e  e 

using  the  same  approach  as  in  subsection  4.2.  An  error  estimate  is  then 

derived.  An  outline  on  how  to  solve  for  yO  using  the  weak  limits  of  the 

e 

eigenvectors  of  A  is  also  given.  •  Let  y^  ■  (y®,y®  )  be  defined  by 
e  e  0  el 


tv 

m 


32  0 
3  y0 


\ 


+  A  yO  =  f  on  Q 


3t2  000 


y®  =  0  on  Z 
0  0 


3y  r 


8v 


0  on  R 


AO 


y0(0)  *  h  on  fl 
0  0  0 


a  0 

?yo 

-  (0)  =  g  on  fi 

3t  0  0 


a2  0 

3  y  « 

_fi  +  tA  vo 

3t2  1  'cl 


\ 


f  on  Q 
1  1 


el 


0  on  Z 


1 


y0  -  y0  on  R  , 


el 


ay 


- —  (0)  »  g  in  n 
3t  1  1 


(5-19) 


(5-20) 


Problem  (5-20)  is  a  nonhomogeneous  boundary  value  problem.  Using  transposi¬ 


tion,  y0  is  defined  in  a  weak  sense  as  in  subsection  4.2. 
si 

In  order  to  derive  an  asymptotic  error  estimate  between  y  and  y^,  they 

e  e 


must  be  redefined  using  transposition.  For  this  purpose,  let 
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.  3*o  3  Xo 

}  xO:  X,.  e  L2(0,T;  H^Cl)),  —  c  l2(Qq)»  — ~  +  AqXq  e  L2(Qq),  Xq  =  0  °n  E^, 


0  /  0 

3X 


3t 


3v 


0  3*0 

=  0  on  R,  x  (x,T)  =  0,  -  (x,T)  «  0 


>1 


Ao 


0 


3t 


*X, 


(-1  d  Xj 

*  *  |x  '•  X  e  L2(0,T;  H^ft  )),  -  e  L2(Q  ),  -  +  e  A  x  e  L2(Q  ),  x  *  0 

1  /  1  1  1  3t  1  3t2  l  l  11 


3X 


1 


1 


on  E  U  R  ,  X  (*,T)  “0,  -  (x,T)  -  0 

1  1  3t 


It  can  be  easily  verified  that  y  ,  yO  satisfy 

e  e 


32X0  \  3X0 

I  y  (  — —  +  A  x  JdQn  ■  /  f  X  dQ  -  /  h  — -  (x,0)d« 
n  e  \3t2  0  0/  0  n  0  0  0  n  0  3t  0 

^0  0 


% 

3y 


el 


+  /  g«X  <x.0)dflrt  +  e  /  - -  X  dR  »  for  all  xn  e  * 

.  t  Qq  0  0  0  R  3VA1  0  0  0 


3  X 


1 


3X, 


/  y  I  -  +  t  A  x  I  dQ  *■  /  f  ft  dQ  -  /  h  -  (rt,0)dfl 

n  ci\3t2  ll/il,  1  11  1  1  1  3t  1 


3X 


1 


+  /  S,X ,(*.  0)ja  -  e'/  y  - -  dR  ,  for  all  x,  o  * 

'•SI,  1  1  1  R  E°  SU*1  1  1 


is  endowed  with  the  topology  carried  over  by  the  mapping 
♦  t  e  L2(0,T;  L2^))  | — ►  Xt ,  i  “  0,1  . 
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(5-21) 


(5-22) 


i 


Kil 


1 


I*2*0  \ 

'2 \TJ  *  Vo)  dQ0  ■ 


=  /  f  x  «*Q  -  /  h  -  (x,0)dQ 

1  0  0  0  :  0  3t  0 

70  0 


(5-23) 


+  /  g  X  (x,0)dfl  ,  for  all  x  e  * 

0  0  0  0  0 

no 


!*\  \ 

J  y0  I -  +  e  A  x  dQ  -  f  f  X  dQ  “  /  h  -  (x,0)d0 

0.  el  W  1  l/  1  o.  1  1  1  o  1  3t  1 


(5-24) 


+  J  g  X  (x,0)dft 
0  1  1  '  1 


r 

e  J  y  -  dR  ,  for  all  x  e  ® 

R  0  3VA!  1  1 


Theorem  5.2: 


Let  y  ,  yO  be  the  solutions  of  (5-21 )— (5— 24 ) .  Then  for  sufficiently 
e  e 


small  e,  the  following  estimate  holds 


ly  -  y0*  <  C  e1/7* 

G  £  L2(Q) 


(5-25) 


Proof: 


Subtract  (5-22)  from  (5-21)  and  (5-24)  from  (5-23)  to  get: 


,  /  0\  /5  X0  \  .  Syel 

J  ( y  -  I  “7“  +  ArtXn  dQn  “  "  £  /  7 -  Xn  dR 

Q0  \  °)  \  3t  00  /  0  R  3Va1  0 


(5-26) 


'i(v-^)(7T*£V,)dV- 


e  /  (y  —  y0) - dR 

R  eo  y0  3vAl 


(5-27) 


+  A  x  -  y  ■  y®  on  Q 
3t2  0  eO  eO  0  0 


Y  =  0  on  E 
eO  0 


=  0  on  R 


X  (x,T)  *  0  on  8 
eO  0 


(x,T)  ■  0  on  ft 
3t  0 


(5-28) 


32X 


el 


— 7*  +  E  A  X  y  "  y°  on  Q 
3t2  1  el  el  el  1 


X  *  0  on  I 
el  l 


X  «  0  on  l 
el 


X  (x,T)  «  0  on  fl 
el  1 


3X 


el 

-  (x,T)  «  0  on  8 

3t  l 


(5-29) 
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It  can  be  shown  (4]  that 


Consequently 


X  e  ♦  .  I  *  0.1  • 

el  i 


eO 


3X 


el 


d\)| 


e  L2(0,T;  H1/2(S)) 


e  L2(0,T;H_1/2(S))  . 


From  (5-13),  one  concludes  that 


e  I  y 


eO 


I  <  C 

R  L2(0,T;  Hl/2(S))  1 


(5-30) 


rri 


3yel 

3vAi 


<  c 


r  L2(0,T;  H"1/2(S))  2 


(5-31) 


Let  x  ■  X  in  (5-26)  and  x  “  X  in  (5-27)  and  use  (5-28)-(5-29) ,  (5-30)- 
0  eO  1  el 

(5-31),  to  get 

ly  -  y0|  <  C  e1/* 

C0  0  L2(Q0)  3 


It  -  *0  I  <  c  e1/4 

“  El  L*(Qj)  4 


from  which  one  obtains  (5-25). 


Now  the  weak  limits  of  the  eigenvectors  of  A  ,  i.e.,  (4-25)-(4-26)  are 

e 

employed  to  solve  (5—1 3)— (5—20) .  First,  renormalize  .  Then,  the 

k-1 

solution  of  (5-19)-(5-20)  may  be  represented  by: 
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twnwifiwiiil  iVf 


yO  »  E  ck  4«k 
0  k=l  0 

®  T 

ck  e  L2(0,T),  E  J  |ck(t)|2  dt  <  » 
k=l  0 


yO  =  i  dk  $k  . 

el  ic*i  e  1 

dk  e  L2(0,T),  E  /  |dk(t)|2  dt  <  -  . 

e  k-1  0  e 


In  order  to  obtain  ck(t),  dk(t),  let 


d0 

X  (x,t)  -  0(tHk(x),  0  e  C2([0,T]),  0(T)  -  ~  (T)  -  0 
0  0  dt 


dv 

X  (x,t)  -  v(tHk(x),  v  e  C2([0,T]),  v(T)  «  7*  (T)  «  0 
1  1  dt 


in  (5-23)-(5-29)  to  get 


d20  \  T  d0 

-  +  yk  0  dt  *  J  (f  ,+k)0dt  -  (h  ,\|»k)  —  (0)  +  (g  ,*k>0(O) 

dt2  0  /  0  .0  0  0  0  dt  0  0 


T  /d2v  \  T  ,  dv 

/  dk -  +  Xk  ev  dt  -  /  (f  ,$k)0dt  -  (h  ,*k)  7-  (0) 

0  e  \dt2  1  /  0  1  1  1  1  dt 


T  /  -  Hi  \ 

+  (g  ,*k)v(0)  -cf  l  C*  /  ♦*  - -  df>  I  vdt 

ll  0  \£-l  S  °  3vAi  / 


B-39 


►M  -Vr 


which  are  equivalent  to 


d2ck 

-  +  pk  ck  =  (f  ,\)ik) 

dt2  0  U  0 


ck(0)  =  (h  ,*k> 
0  0 


dck  . 

7—  (0)  =■  (g  ,*k) 
dt  0  0 


-  +  Ak  e.  dk  -  (f  ,$k) 

dt2  1  e  11 


-  3*1 

e  £  c^  J  - -  dS 

1-1  S  0  3\>A! 


dk(0)  -  (h  ,*k) 
e  1  1 


d  d* 

(0)  -  <g.*k) 
dt  11 


Remark  5.1: 

It  should  be  stressed  once  more  that  the  concepts  used  in  approximating 
( 5—5 )— (5— 1 1 )  can  be  employed  to  approximate  any  boundary  value  problem  that 
fits  into  (5-l)-(5-5). 
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SECTION  6 


CONCLUSIONS 

In  this  paper,  three  stiff  classical  boundary  value  problems  (namely 
elliptic,  parabolic  and  hyperbolic)  have  been  analyzed,  using  the  spectral 
analysis  of  stiff  operators  [5]. 

For  elliptic  problems,  it  has  been  shown,  that  by  appropriately  modifying 
the  weak  limits  of  the  eigenvectors  of  Ae,  a  Laurent  series  expansion  of  any 
order  can  be  derived.  For  evolution  problems,  zeroth  order  approximations 
in  L2(0,T;  L2(0))  were  easily  constructed,  using  only  the  weak  limits  of  the 
eigenvectors  of  Ag. 

Many  other  similar  problems  (including  control  problems)  can  be  analyzed 
using  the  concepts  used  in  the  present  paper. 

Nice  physical  interpretations  can  be  associated  with  the  behavior  of  the 
solutions  of  (l-l)-(l-3)  as  e  +  0  in  the  fields  of  heat  transfer  and  electro¬ 
magnetic  wave  propagation.  In  a  future  publication  [6],  these  interpretations 
are  discussed  using  several  examples  of  stiff  boundary  value  problems. 

The  approximations  obtained  in  this  paper  are  also  compared  with  finite 


dimensional  approximations 
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APPENDIX  C 


C.l  INTRODUCTION 

four  different  models  have  been  used  during  the  design  and  analysis  of 
the  report*  The  models  are  the  subsystem  design  model,  the  subsystem  truth 
model,  the  coupled  system  analysis  model,  and  the  coupled  system  truth  model. 
Each  model  is  a  low-order  finite  element  approximation  of  a  set  of  masses 
corrected  by  flexible  beams.  The  models  take  the  form: 

Mp  +  Dp  +  Kp  -  Bf  (C-l) 

y  -  Cp  (C-2) 

where  N  ■  number  of  elements 
m  -  number  of  forces 

p  -  N  vector  of  positions  of  the  elements 
f  -  N  vector  of  forces  or  the  elements 
y  *  r  vector  of  measurements 
M  *  NxN  mass  matrix 

» 

K  -  NxN  spring  constant  matrix 
K  ■  NxN  damping  matrix 
B  -  Nxm  force  coefficient  matrix 
C  ■  rxN  output  matrix 

The  remaining  sections  present  the  data  for  each  of  the  finite  element  models 
that  has  been  used. 
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C.2  SUBSYSTEM  DESIGN  MODEL 
Number  of  masses  *  3 
N  =»  5 


M  - 

5.5000 

0.3536 

0.0000 

0.3536 

5.5000 

0.3536 

0.0000 

0.3536 

5.5000 

K  - 

5.3924 

-3.8130 

0.0000 

-3.8130 

5.3924 

-3.8130 

0.0000 

-3.8130 

5.3924 

D  * 

0.3300 

0.0212 

0.0000 

0.0212 

0.3300 

0.0212 

0.0000 

0.0212 

0.3300 

B  - 

2.1213 

0.0000 

0.0000 

0.0000 

2.1213 

0.0000 

0.0000 

0.0000 

2.1213 

C  - 

1.0000 

0.0000 

-0.7071 

0.7071 

0.0000 

-1.0000 

Note:  The  inputs  for  subsystems  1  and  2  are  defined  in  Fig.  2-1. 

C.3  SUBSYSTEM  TRUTH  MODEL 
Number  of  massess  *  3 
N  -  5 


M  - 

B.5000 

0.3536 

0.0000 

0.0000 

0.0000 

0.3536 

1.0000 

0.2500 

0.0000 

0.0000 

0.0000 

0.2500 

B.5000 

0.2500 

0.0000 

0.0000 

0.0000 

0.2500 

1.0000 

0.3536 

0.0000 

0.0000 

0.0000 

0.3536 

B.5000 

K  - 

18.7500 

-13.2583 

O.'OOOO 

0.0000 

0.0000 

-13.2583 

18.7500 

-9.3750 

0.0000 

0.0000 

0.0000 

-9.3750 

18.7500 

-9.3750 

0.0000 

0.0000 

0.0000 

-9.3750 

18.7500 

-13.2583 

0.0000 

0.0000 

0.0000 

-13.2583 

18.7500 

D  - 

0.5100 

0.0212 

0.0000 

0.0000 

0.0000 

0.0212 

0.0600 

0.0150 

0.0000 

0.0000 

0.0000 

0.0150 

0.5100 

0.0150 

0.0000 

0.0000 

0.0000 

0.0150 

0.0600 

0.0212 

0.0000 

0.0000 

0.0000 

0.0212 

0.5100 

C-2 


0,0000 

0.0000 

0.0000 

0.0000 

2.7386 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.7386 

C  =  1.0000 

0.0000 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.7071 

0.0000 

-1.0000 

COUPLED  ANALYSIS  MODEL 

Number  of  masses  *  6 

N  -  6 

M  *  5.5000 

0.3536 

0.0000 

0.0000 

0.0000 

0.0000 

0.3536 

5.5000 

0.3518 

0.0000 

0.0000 

0.0000 

0.0000 

0.3518 

5.5000 

0.0025 

0.0000 

0.0000 

0.0000 

0.0000 

0.0025 

5.5000 

0.3518 

0.0000 

0.0000 

0.0000 

0.0000 

0.3518 

5.5000 

0.3536 

0.0000 

0.0000 

0.0000 

0.0000 

0.3536 

5.5000 

K  -  5.3924 

-3.8130 

0.0000 

0.0000 

0.0000 

0.0000 

-3.8130 

5.3924 

-3.7941 

0.0000 

0.0000 

0.0000 

0.0000 

-3.7941 

5.3924 

-0.0270 

0.0000 

0.0000 

0.0000 

0.0000 

-0.0270 

5.3924 

-3.7941 

0.0000 

0.0000 

0.0000 

0.0000 

-3.7941 

5.3924 

-3.8130 

0.0000 

0.0000 

0.0000 

0.0000 

-3.8130 

5.3924 

D  -  0.3300 

0.0212 

0.0000 

0.0000 

0.0000 

0.0000 

0.0212 

0.3300 

0.0211 

0.0000 

0.0000 

0.0000 

0.0000 

0.0211 

0.3300 

0.0002 

0.0000 

0.0000 

0.0000 

0.0000 

0.0002 

0.3300 

0.0211 

0.0000 

0.0000 

0.0000 

0.0000 

0.0211 

0.3300 

0.0212 

0.0000 

0.0000 

0.0000 

0.0000 

0.0212 

0.3300 

B  -  2.1213 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.1213 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2*1213 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.1213 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.1213 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.1213 

C  -  1.0000 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.7071 

-0.9950 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.9950 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-0.9950 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.9950 

0.7071 

-1.0000 
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C. 5  COUPLED  TRUTH  MODEL 


Number  of  masses  =  6 


N  =  11 


M 


K 


D 


B 


C 


9.2500 

0.3536 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.3536 

1 .0000 

0.2500 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.2500 

9.2500 

0.2500 

0.0000 

0.0000 

0.0000 

0.0000 

0.0Q00 

0.0000 

0.0000 

0.0000 

0.0000 

0.2500 

1 .0000 

0.3518 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.3518 

9.2500 

0.0025 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0025 

1 .0000 

0.0025 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0025 

9.2500 

0.3518 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.3518 

1.0000 

0.2500 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.2500 

9.2500 

0.2500 

O.ftOOO 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.2500 

1.0000 

0.3536 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.3536 

9.2500 

18.1244 

-12.8159 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

' 

0.0000 

-12.8159 

18.1244 

-9.0622 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

O.OOOC 

0.0000 

-9.0622 

18.1244 

-9.0622 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-9.0622 

18.1244 

-12.7523 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-12.7523 

18.1244 

-0.0906 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-9.0906 

18.1244 

-0.0906 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-12.7523 

18.1244 

-12.7523 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

1  0.0000 

0.0000 

0.0000 

0.0000 

-12.7523 

18.1244 

-9.0622 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-9.0622 

18.1244 

-9.0622 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-9.0622 

18.1244 

-12.8159 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-12.8159 

18.1244 

0.5550 

0.0212 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0212 

0.0600 

0.0150 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0150 

0.5550 

0.0150 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0150 

0.0600 

0.0211 

0.0000 

0.0000 

0.0000 

0.0000 

o.ocoo 

0.0000 

0.0000 

0.0000 

0.0000 

0.0211 

0.5550 

0.0002 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0002 

0.0600 

0.0002 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0002 

0.5550 

0.0211 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0211 

0.0600 

0.0150 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0150 

0.5550 

0.0150 

0.0000 

0.0000 

o.ooco 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

o.oiso 

0.0600 

0.0212 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0212 

0.5550 

2.8/23 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.8723 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.8723 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.8723 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.8723 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

> 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

2.8723 

1.0000 

0.0000 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.7071 

0.0000 

-0.9950 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.9950 

0.0000 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-0.9950 

0.0000 

-0.7071 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

0.9950 

0.0000 

0.7071 

0.0000 

-1.0000 
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